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Ne Om ene 
A. DIGITAL SIGNAL PROCESSING 


Digital signal processing has ie HOOtG Ih, UGthe@enituny, 
mathematics especially in the fields of astronomy arfd com- 
pilation of mathematical tables. Today it has become a 
powerful tool in a multitude of diverse fields of science 
and technology. The application of digital signal process- 
ing varies from low-frequency spectrum seismology through 
spectral analysis of speech and sonar into the video spectrum 
of radar systems [1]. 

Rabiner and Gold have interrelated in their excellent 
book on "Theory and Application of Digital Signal Processing" 
(see Ref [2]) digital processing theory, 
with a variety of applications ranging from radar, sonar, 
communication, music, seismic and medical signal processing, 
and, with digital component technology. This technol- 
ogy is the main driving force for progress in this field 
as well as the general area of computer design. It is also 
observed by them that although the formulation of engineer- 
ing problems is often as vague as those of the "softer" 
Ssclences (such as anthropology, psychology, etc.) the execu- 
tion of these problems appears to depend on greater and 
greater accuracy and reproducibility. The capability of 
digital systems to achieve a guaranteed accuracy and essen- 


tially perfect reproducibility is very appealing to engineers. 


ie 





A digital filter is defined in [3] to be a computational 
process or algorithm by which a digital signal or equiva- 
lently a sequence of numbers acting as input, is transformed 
into a second sequence of numbers or output digital signal, 
where. the term 'digital' implies that both time (the inde- 
pendent variable) and amplitude are quantized. A digital 
Signal 1s said to be a discrete signal whenever the ampli- 
Eude quantization step of a time sampled signal is zero. 

The field of one-dimensional digital filtering, which is 
outlined in Figure 1.1, encompasses recursive, nonrecursive 
and Fast-Fourier Transform processing. The terms recursive 
Bagoemenrecursive, instead of IIR (Infinite Impulse Response) 
comes FIR (Finate Impulse Response) are preferred in the con- 
text of this realization oriented thesis. Kaiser showed in 
[4] that recursive processing is much more efficient than 
nonrecursive processing. Stockham's [5] method to perform 
fast convolution, which later became known as FFT method, 
improved the efficiency of nonrecursive techniques, so that 
comparisons in one-dimension are no longer strongly biased 
toward recursive methods [2]. 

One-dimensional recursive digital filter design depends 
stongly on the effective and well developed continuous fil- 
ter design theory. Moreover, the stability analysis of 
higher order one-dimensional recursive realizations can be 
solved by investigation of the root distribution of its fac- 
tored form representation. It is important to realize that 


continuous domain design techniques and the factorization 


14 





property (the fundamental theorem of algebra) exist only in 
one-dimension. 

Areas of applications of one-dimensional digital process- 
ing are listed in Figure 1.1. There are important military, 
and specifically naval applications in addition to those 
listed, such as missile and torpedo guidance, launcher con- 
trol system, fire control, combat information collections, 
dissemination and display. 

Digital filtering in several dimension, which 1s over- 
viewed in Figure 1.2, has gained considerable importance, 
especially for the two-dimensional case. Figure 1.3 displays 
an example of two-dimensional low-pass and high-pass filter- 
ing. Until 1966, two-dimensional digital filtering was imple- 
mented by nonrecursive or convolutional techniques. In this 
method the output is the weighted sum of unit sample responses 
of all past input values. The serious disadvantage of the 
convolutional method is the requirement of a very large num- 
ber of arithmetic operations. The development of the Fast- 
Fourier-Transform (FFT) algorithm in 1966, reduced the number 
of arithmetic operations considerably and is used extensively 
today. Filtering via FFT is accomplished by computing the 
transform of the input function, multiplying by the filter 
transfer function, and inverse transforming the result. The 
recursive algorithm, which has in general an infinite impulse 
response [3], constitutes another technique for the realiza- 
iron Of two-dimensional digital filters. 

Hall [6] compares the amount of computation required for 


the three filtering techniques. He shows that, in general, 
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FIG 1.3a EXAMPLE OF TWO-DIMENSIONAL DIGITAL FILTERING, [2] 
CORIGINAL PHOTOGRAPH ) 
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FIG 1.3b EXAMPLE OF TWO-DIMENSIONAL DIGITAL FILTERING, [2] 


(LOW-PASS FILTERED VERSION OF ORIGINAL PHOTOGRAPH ) 
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PHOTOGRAPH ) 


AO 





the FFT and recursive algorithm are preferable to the non- 
recursive oneinnumber of computations and storage ce 
ments. Hall demonstrates that the recursive filter algorithm 
constitutes the best method for large amounts of data, l.e., 
is the fastest and the cheapest. 

When processing two-dimensional data by recursive filters 
a fundamental problem exists due to the inherent feedback, 
namely, the problem of numerical stability. Since, in 
general, two variable polynomials eee be factored into a 
product of real first and second order real coefficient poly- 
femials;, it 1s difficult to solve the stability problem for 
two-dimensional recursive filtering. Consequently, the 
majority of papers published in two-dimensional digital fil- 
tering deal with the design of nonrecursive filters which 
are inherently stable, i.e., [7] - [11]. 

There are several reports discussing two-dimensional 
Peeursive digital filters [12]-{15]. For example Miii@and 
[15] formulate two-dimensional recursive filters by the Z- 
transform and linear difference equations, and although they 
investigate problem areas related to stability and realiza- 
tion, the majority of problems remain to be solved. 

The most important applications of two-dimensional digi- 
tal filtering of digital data are in picture processing and 
geophysical data analysis. Picture processing can be cate- 
Porazed into; 

ales Deedee image restoration and enhancement , and 


2. Computer pictorial pattern recognition. 
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The methods of image restoration and enhancement are 
applied to invert degradations, such as aberrations, atmos- 
pheric effects, scanning, motion, and to manipulate images 
to improve viewing phenomena experienced by the human eye. 
Image restoration and enhancement has been used with great 
success in biomedicine [16], i.e., in extraction of quanti- 
tative information from X-ray films, chromosome counting, 
measuring the extent of atherosclerosis from arteriograms, 
in forensic sciences application, i.e., fingerprint image 
enhancement for automatic classification, and in astronomy 
[17], i.e., removal of turbulence from astronomical photo- 
graphy [18]. Pictorial pattern recognition by digital com- 
puter has gained great importance in satellite surveillance 
and military reconaissance. The application of two-dimen- 
Sional digital filter permit the separation of different 
horizontal scans on magnetic and gravity maps in geophysics 
and can be used equally well for structural and topographic 
maps or for any other type of data which is available in 
the format of a planar grid [19], [4]. 

Although digital filtering in several dimensions is sane 
ing importance in medicine, i.e.,heart volume measurements, 
in fire control problems, and to analyze complex electronic 
circuitry, there exists no general theory concerning struc- 
tures, analysis and design. 

These considerations are summarized as follows: since 
one-dimensional recursive digital analysis and design tech- 


niques depend strongly on the fundamental theorem of algebra 


ae 





and on extensive utilization of continuous design theory, it 
represents a special, i.e., a non-generalizable field in the 
area of N-dimensional recursive digital filtering. Multi- 
dimensional digital filtering is well developed in two-dimen- 
Sions but is based, due to the absence of recursive theory 
analysis and design tools, predominantly on nonrecursive 


filtering schemes. 
B. AREA OF INVESTIGATION 


The design of N-dimensional recursive digital filters is 
achieved by choosing filter coefficients to approximate given 
Specifications and to evaluate the stability of the resulting 
filter. Several methods, which are discussed in Chapter IT, 
are available to determine whether a multi-dimensional recur- 
Sive digital filter is stable. They are based on a numerical 
root distribution investigation of the N-dimensional charac- 
teristic equation of H(z) and determine whether or not a fil- 
ter is stable.” There is no qualitative information extracted 
about how stable or unstable a system is or what coefficient 
changes of H(z) are necessary to make an unstable system 
stable. Moreover, these methods are complicated and result 
in time-consuming computer programs, which must be applied 
to each individual filter design variation. Due to these 
reasons, the design of multi-dimensional, recursive digital 
filters has not been achieved [2]. 

In the following chapters, a different approach to anal- 
ysis and design of N-dimensional recursive filters is pre- 


sented which is based on the Taylor expansion of the filter 
The notation is discussed in Appendix A. 
2s 





H(z). It can be shown that an N-dimensional recursive filter 
1s stable if an only if the sum of all unit sample response 


entries converge absolutely, i.e., 


S = y petGa, ) || ee ee 


ni 

(O<n,<@) 
where the h(n) are the coefficients of the Taylor expansion 
of N-dimensional transfer function H(z). 

Thus, conditions for stability of an N-dimensional recur- 
sive system H(z) can be obtained by observing the absolute 
convergency behavior of its Taylor expansion coefficients. 
Although this represents a rather complicated method to find 
stability conditions, it can be used to identify important 
unit sample response characteristics of certain two-dimen- 
Slonal recursive filters, which in turn can be generalized 
to N-dimensions. 

The relationship of the N-dimensional recursive filter 
coefficients to the unit sample response coefficients is 
investigated and methods to perform the approximation step 
of time-domain design are studied. Also, structural realiza- 
tion of N-dimensional transfer functions are searched. 

In summary a new design and implementation technique for 


N-dimensional recursive filters is developed. 
fer ke VIEW OF RESULTS 


The results can be grouped into techniques for stability 


analysis of N-dimensional recursive filters (Chapters V to 
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VII), three techniques to solve the approximation part of 
N-dimensional recursive filter design (Chapters VIII and IX) 
and, combining both previous results, design techniques are 
derived (Chapter X). In addition, the implementation of 
N-dimensional recursive filters is considered and some new 
forms presented. These results are obtained using a compact 
vector index notation and a general theory defining charac- 
teristics of N-dimensional systems, as well as the N- 
dimensional Z-transform and some of its properties (see 
Appendix A). 

The computations of N-dimensional Taylor expansion coef- 
fierentcen(n) are discussed in Chapter IV. The standard 
derivative operator method requires coneneieteed computations, 
and the convolutional or recursive technique requires the 
knowledge of all previous values to compute h(n). Methods 
of combinatorial algebra are applied to derive two efficient 
nonrecursive formulas for Taylor expansion coefficients, of 
which one has high speed characteristics, 1.e., requires very 
few multiplications and additions. The results of Chapter IV 
are used in Chapter V to solve the stability problem of the 
second-order (coupled), third- and fourth-order (uncoupled) 
two-dimensional recursive digital filter using Ultraspherical 
and derivative operator theory. 

In Chapters VI and VII, observations made in two- 
dimensional stability analysis of the previous chapter are 


generalized to cases of any order in two-and also in 
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L—Gemensions, resulting in new sets of sufficient and 
necessary stability conditions. <A general method to. derive 
necessary and sufficient stability conditions for N-dimensional 
recursive digital filters is presented and new stability 
conditions for the first degree and positive coefficient case 
in N-dimensions are listed. 

The important advantage of the new approach to the 
stability problem over the previously outlined numerical ones 
becomes evident if one considers that for a given type of 
filter transfer function the methods of Chapters V to VIl are 
eamxoried Only once to obtain stability @eonditions in terms of 
eoefficients of the characteristic equation of H(z). These 
conditions can be used by the design engineer in a simple 
manner to test whether a specific filter is stable. They 
elow tne identification of coefficients causing instability 
pugmean casily be interpreted to what extent coefficients 
must be changed to stabilize an unstable system. 

The nonrecursive combinatorial formula of Chapter IV is 
used in Chapter VIII to formulate rules for the propagation 
of the characteristic equation coefficients in H(z) in the 
Gerresponding unit sample response. Conversely, computer 
algorithms for the extraction of an exact all-pole or 
meetonal filter transfer function or a parallel arrangement 
of low-order all-pole or rational filter from a given unit 
Sample response are presented in Chapter IX. All of these 


methods are new. The extraction of a parallel arrangement 
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of low-order sections can be viewed as an approximate method 
for partial fraction expansion in N-variables. 

In Chapter X, a time-domain design technique is presented 
using an extraction technique to approximate a given specifi- 
cation by a transfer function. The stability of the extracted 
transfer function is analyzed using cava Hercule presented in 
Chapters V to VII. A measure of quality for the N-dimensional 
eerroxlmMacaon techniques 1s introduced which compares the 
spectra of the given unit sample response with the spectra 
@meuie Extracted transfer function. Finally, the implications 
of the N-variable partial fraction expansion technique on 
frequency domain design are discussed and areas of fruture 
mesearch are outlined. 

in summary a series method to analyze N-dimensional 
recursive digital filters is presented and new results 
regarding the stability of N-dimensional recursive digital 
filters are achieved. The series approach is used to formu- 
mee tmree new transfer function extraction techniques which, 
combined with the stability and structural analysis, present 
anew time-domain design technique for N-dimensional 


recursive filters. 
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eo eenOUS LO enter SABIE Y OF 
MULPTI-DIMENSTONAL SMSTEMS 


A. INTRODUCTION 
Steiglitz [21] has demonstrated by the construction of a 


specific isomorphism between the Hilbert Space of a complex 
valued Lebesque measurable function f(t) and the Hilbert 


Oo 


Space of double-ended sequences of complex numbers eo eer: 
that the theories of one-diemnsional signal processing with 
linear time-invariant realizable filters are identical in 

the continuous and discrete case. It is interesting to note 
that optimization problems involving quadratic cost functions 
are also equivalent in the continuous and discrete time do- 
mains. 

The extension of this theory to several dimensions, 
though mathematically feasible, has little importance since 
a powerful continuous filter theory exists only in one di- 
mension. 

Another important property of one-dimensional theory is 
deduced directly from the fundamental theorem of algebra 
which states that every higher-order polynomial can be fac- 
tored into a product of lower-order polynomials. The impor- 
tance of the one-variable factorization property cannot be 
minimized. It permits the designer to investigate stability 
of a higher-order system by determining the factored form 


representation and applying simple stability criteria to 


each factor. If each and every factor is stable, then the 
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higher-order system is stable [1], [22]. The factorization 
property also constitutes a method of realization of higher- 
order filters by paralleling, or cascading lower order sec- 
tions, thereby reducing coefficient seeds eA tsp problems 
inherent in realizing higher order systems. Finally, an un- 
stable one-variable filter can be stabilized by cascading 

an appropriate all-pass filter [2]. This does not work in 
tne analog domain because exact pole-zero cancellation is 
impossible. This technique would seem to be feasible in the 
discrete time domain except that roundoff effects may cause 
difficulties. It is shown in Appendix B that multivariable 
polynomials cannot be factored in general. The above are 
some of the major difficulties opposing the generalization 
of one-dimensional system theory to multi-dimensions. Some 
important conceptual ideas of one-dimensional analysis that 
can be generalized to several dimensions are presented in 
Appendix A. These include the definition of N-dimensional 
linear, shift-invariant, realizable, and separable systems, 
the N-dimensional Z-Transform and its properties, as well 

as the notion of the N-dimensional recursive equations and 
transfer RENO r, 

The following paragraphs deal with existing methods to 
determine stability of N-dimensional recursive systems, H(z) 
(notation defined in Appendix A) which are obtained by apply- 
ing the N-dimensional Z-Transform to a linear, shift-invari- 


ant, causal N-dimensional difference equation. 
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Pee vEenobs. LOO DETERMINE SPABILITY OF MULTI-VAREABLE SYSTEMS 


Shanks, Treital and Justice recently pulished a tech- 
nique to determine the stability of a two-dimensional recur- 
PPeveeSyotems H\z,52,), by which the unit disk in the 2)- 
plane must be mapped into the z,-plane by solving the impli- 
cit two-variable denominator polynomial of GZ, 2 9) nec 
essary and sufficient condition for stability is determined 
if the Pitot the Zz splkame unit disk an the z>=plane does 
not overlap the unit disk in the z,-plane [12]. This method 
requires an infinite number of mappings and cannot, there- 
fore, be applied exactly. Huang [13] simplified the stabil- 
ity problem by the observation that 


A(€z15Z2) 


MSDS Garam 


where A and B are two-variable polynomials, is stableif and 
ommey if , 


1. The map of 6d; ={z,:|2,]=1} in the z2-plane accord- 


Mier to the equation B(z,,z.) = 0 lies inside d> = {z2:|z2|21!. 
VeeeNe point im ad, = {z,:|z,| 21) maps imto the point 
fp 0 by the relation BCz,,2,) = 0. 


Although this procedure is finite, it requires the applica- 
tion of two bilinear transformations to formulate the prob- 
lem as discussed by Ansell [23] whose main contribution is 
to couple the use of a Hermite test for checking stability 
with a series of Sturm tests checking positivity. Anderson 
and Jury [24] modified Huang's method by outlining a pro- 


cedure which replaces the bilinear transformations by the 
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construction of a Schur-Cohn matrix and the checking for 
positivity of a set of self-inversive polynomials. It also 
replaces the Hermite test [25] by a Schur-Cohn matrix test 
and requires a series of Sturm tests or equivalently, a 
system of tests establishing the root distribution of a 
polynomial. These modifications represent a substantial 
reduction in computations as compared to Huang's method. 
Siljak [26] comments on [24] by proposing a substantial 
Simplification of the number of Sturm tests which he states 
1s desirable because a computer processing of the Sturm test 
1s not available. However, a suitable computer algorithm 
for the Sturm test has recently been published by Singh and 
Panwar, {27]. The determination of stability of multivari- 
able filters is solved by Bose and Kamat [28] using a method 
outlined in [29] where the N-variable characteristic equa- 
ftemeot H(Z) is reduced to several single-variable polynomi- 
als by an algorithm composed of a finite number of rational 
operations, followed by repeated factorization of the set 
of one-variable polynomials and back substitutions. This 
procedure 1S very complicated and the implementation on a 
general purpose digital computer may not be feasible [28], 
as recognized by the authors themselves. A slightly simpler 
method for determining the stability of two-dimensional 
filters has been offered very recently by Silja Oe ein an 
algorithmic approach requiring the equivalence of two Marden 


tests [29] and checking positive definiteness with methods 
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developed in [31], [32], and [33]. This last approach has 
been implemented by Schaldach [34] on a general purpose com- 
puter. Finally, Jury [36] presented an excellent summary 
on methods to determine stability of multivariable filter 


uSing the inners approach. 


ee SUMMARY 


The above outlined methods to determine the stability of 
multivariable characteristic equation of H(z) to identify 


roots not satisfying 


() 2a gaal C2) 


a= 7): 


were CY Jacl <2 = (Ia) 2)M (la) <1)M--- Mle, I< 2) 


They also have in common the property that they are highly 
theoretical and very difficult to apply [2] and of limited 
practical value to the design engineer. It is necessary for 
analysis and design purposes of multivariable systems to 
identify not only whether or not a system is stable, l.e., 
that all roots of the N-variable characteristic equation 
satisfy eq (2.1), but also what coefficients of an unstable 
transfer function cause the instability and how the coef- 
ficients must be changed to achieve stability. This latter 
question can only be solved definitely utilizing previously 
considered methods for the bilinear two-dimensional case, 


for example: 


Be 
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(2.2) 
HoOrmeWiich Stability Conditions are derived in [12], [13], 
and [24]. This result will be discussed further in Chapters 


weeamnd VIT. 
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IIL. STRUCTURES AND STATE SPACE FORMULATION 
FOR N-DIMENSIONAL DIGITAL FILTER 


A. INTRODUCTION 

There are several ways of realizing a given digital filter 
transfer function. The choice of structure depends on finite 
wordlength effects in practical hardware realizations on 
eemplexity of computation. 

Sie-dimensional canonical structures are most economical 
and fastest, yet they may have undesirable finite register 
length effects. In this chapter a one-dimensional direct 
form filter realization is introduced designated by direct 
form 4 to distinguish it from other direct forms which 
corresponds directly to hardware implementation and by 
realizing the numerator and denominator polynomial separately 
in cascaded form. The direct form 4 is deductively used to 
derive direct form and canonical structures in N-dimensions, 
which have considerable importance, since the nonfactora- 
bility characteristic of multivariable polynomials eliminates 
the possibility of realizing high-order multi-dimensional 
filter by paralleling or cascading low-order standardized 
sections, except for the separable case. 

The state space approach provides a means for character- 
izing complex systems in a simple manner. The complexity of 
manipulation is decreased and the easy implementation on a 


general purpose digital computer minimizes computational 
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difficulties. The increase of the number of state variables; 
maieenumber Of Inputs, or the number of outputs does not in- 
crease the complexity of the equations. In fact, the analysis 
Got complicated multiple-input-multiple-output Systems can be 
eienied OUL by procedures that are only slightly more compli- 
cated than those required for the analysis of systems of first- 
order equations. 

ier Oone-dimensional State variable approach to digital 
filtering was presented by Grindon [37] in 1970. Since then it 
is commonly used in one-dimensional digital filter analysis. 
For example, Parker and Hess [8a] used the state space approach 
Memo@etive twenty four second-order canonic digital filter seéec- 
Eons . 

Parker, Girard and Souchon [42] used state space methods 
Memerepress the co-variance matrix of state errors and the 
meenenes Of the Output Error in terms of the correlation ma- 
memeese Or the Quantization Error sources. This method allowed 
m@@emecvaiUation of correlation effects between error sources in 
Pmemeonparison Of different filter structures. 

Parker and Yakowitz [44] used the state space approach to 
compute quantization error bounds. 

It can be summarized that the state space approach is a 
powerful and commonly used method to analyze one-dimensional 
digital systems. Since the following chapters are mainly con- 
cerned with N-dimensional filter theory, an extension of the 
One-dimensional state space analysis to N-dimensions is devel- ~ 


mead in Section C of this chapter. 
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bee oT RUCTURES OF N-DIMENSZONAL DIGETAL FILTER 


1. One-Dimensional Digital Filter Realization 


The rational transfer function (A.10) has for the 


eme—dimensional case, 1.e., N = 1, the form: 


L 

wy By ra 

ss a 

; = ; An) 

FIC Zae) SS a Se on ea ea CAD) (3.1) 


—) 


! 
= 
=" 
N 
a. 


The time domain representation is obtained by applying the 


inverse Z-Transform to eq (3.1) which results in: 


M 


Ib: 
yien, ) = » oa SAC TR Sn pe B, x(n,-2) C352) 


m= kb L=0 


The realization implementing eq (3.2) directly, is known as 


direct form 1 and is shown in Figure 3.1. 
Equation (3.1) can be rewritten ina slightly different 


form by introducing a new variable, W(z,), such that 


W(z1) ¥Cz)) 
H(z) = X(z,) Wz,) (3.3) 


The corresponding set of difference equations consists of: 
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Winn = x(n es a  w(n,-m) (3.4) 


ae 


= 
I 
}-4 


and 


Bp w(n,-2) Cas) 


L 
ce 
=Q 


) 


Equcrtons (3.4) and (3.5) are mee in the direct form 
2 as shown in Figure 3.2 which can be visualized as a cascade 
realization of two digital filters, realizing the denominator 
and numerator polynomials, respectively. Figure 3.2 can be 
redrawn as shown in Figure 3.3. The resulting realization 
1s called direct form 3 or canonic form, since it represents 
a structure with the number of delays equal to the order. 

At this point, an additional realization will be intro- 
duced, the direct form 4, which combines the characteristics 
of the direct forms 1 and 2 by having only two entries in 
each summer, which corresponds directly to hardware implemen- 
tation and by realizing the numerator and denominator poly- 
nomial separately in cascaded form. The usefulness of the 
direct form 4, which is shown in Figure 3.4 and which re- 
places the direct forms 1 and 2, will become evident in the 
following section. 

There are several other methods to realize a one-dimen- 
sional mth order digital filter. For example, the cascade 


form of first and second order sections H(z,), as shown in 


Sif 





Figure 3.6, where 


k 
Zest Se Gze k = Ea 


Ae (3.6) 


and [x] = largest integr < x 


Each component H. (2, ) can be realized in one of the above 
outlined forms. 

Partial fraction expansion methods applied to eq (3.1) 
lead to a parallel from realization of first and second 


order sections (see Fig 3.5) l.e., 


k 
H(z,) = C + » Hi (21) k = ES} 7) 


1=1 

Other forms of realization include hybrid structures, l.e., 
parallel-cascade forms, transpose configurations, which can 
be obtained for all previously outlined structures by revers- 
ing the direction of signal flow and by interchanging all 
branch and summing modes, wave-digital filter [a6 | to 71500 
and continued fraction expansion filter [51]. 

The direct forms discussed above and the series/parallel 
arrangement of lower order sections are by far the most often 


used ones in computer simulation and in digital hardware. 


2. Deductive Approach to Develop N-Dimensional Structures 


The fundamental theorem of algebra does not extend 
to the multi-variable case which 1S proven in Appendix B. 
Consequently, the series and parallel realization schemes 


cannot be generalized to N-dimensions. The following 
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paragraphs develop a generalization for direct forms to two- 
dimension and then to N-dimensions. 


The numerator of the two-dimensional transfer function 


L, 2 
—X, me 
rh Omens. ae 
2,20 £420 
H(z.,z,) = —_—_-—CS (3.8) 
oe i, 
)2e OL zg, ig, Mme 
m, Mm, l 2 
m, =0 m, =0 
(m,+m,#0) 


can be rewritten as one-variable polynomial in z,, weighted 


by coefficients which are functions in 2,, 1l.e., 


L L 


2 1 
N(z,,Z.) = » » By 9, 21 Teles 2 (3.9) 
; 8-0 .\2,=0 
Ly 
Z > N aye lee 2 4.109 
ps ue oe ? 
25=0 


The same reasoning applied to the denominator polynomial 


leads to 
M, M, 
DCZ,,;22) = lL = ». Yon m2 Za 2 C311, 
, m,=0 m,=0 


_—- =» sce 





Mh 
ie iy (Preinf)) (Bel 


m, =0 


where the sum of the indices is mi + mo # 0. To obtain 
the direct form 3 representation, a new function, W(z,,zZ,), 


mmantroduced: 


W(Z15Z2) n#OZ15 Zo) 
FC Zs Sa) = Gama.) : W(z15Z2) C3 ei30 


where each factor can be written using eq (3.10) as 


Lig 
WHGZ 5 Zine) -) (x CZ) Zz? AT CAR) CS) 
ee) “2 
R5=0 
and eq (3.12) as 
ul 
¥(z,,2,) = W(z,,2,) + yy (Ds 5 2,3 | ap ON Canebapyy (3, 153) 
tlie 
M = 0 


The realization of eqs (3.14) and (3.15) is shown in Figure 
3.7 using compact notation, and in Figure 3.8 where 

Ee, (z,) and Pum 64? are implemented for each (L,£,) and 
and (Mm) by one-dimensional direct form 4 realizations. 
The direct form 3 (canonic) will be derived in the next 
@eetion. The direet form ¥ realization of H(€z,,z,) is 


Shown in Figures 3.10 and 3.11, using compact notation and 


one-dimensional direct form 4 implementation of each Nig 
1¥2 
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and D > respectively. The generalization to N-dimensions 


Mm, 


follows directly where the N-dimensional transfer function 
is written in the numerator and denominator as a one-vari- 
able summation which is weighted by functions of (N-1) vari- 


ables, il.e., 


L 


N 
~hy 
( x, Ne cule Q 21s tyer EN 
i N-1°N 
LW =0 
H(Z1,+++ 52, ) = EE CaaS) 

N 
: = 1y 
1- » (dy a my (219+ +9242] z,, 

m7 0 N-1 


The direct forms 3 and 4 in N-dimensions are shown using com- 


pact notation in Figures 3.12 and 3.13, respectively. 


Example 3.1: 

The direct form 4 realization of a full quadratic trans- 
fer function in three-dimensions is shown in Figure 3.15. 
Structures of two-dimensional third-order and fourth-order 
uncoupled, as well as second-order and fourth-order coupled 
filters, which will be investigated in detail in the next 
chapters, are presented in Appendix G. 


A separable system is defined in Appendix A(4) by 


N 
fags) me | h,(n,) Comey) 


which is shown in Appendix A(5) to imply: 


Ky 





N 
H(z) = M H,(2;) (3.18) 


But H(Z) can be realized by a cascade structure of one-dimen- 
Sional transfer functions in each dimension, as shown in 
Figure 3.16, by techniques discussed in Paragraph 1 of this 


Beet lon. 
SeeeolkATE EQUATION FORMULATION 


1. State equation formulation in two-dimensions 


The two-dimensional transfer function H(z,,z,) can be 


written in long hand notation as 


Ly Ly 
Seats 
2, 2.82, 0 
2, =0 Lo=0 
H(z,,2,) = ——S—————————C—O:-.- 199) 


m=O m,=0 


(m,+m, #0) 


ieee runetion W(z,z,) 1s introduced similarly, as in the 


previous section, then 


W(z,,2,) Y(z,,z,) 
H(z,,Z,) = ———— * Gace 
19“2 WiezA5z >) W Zino» 


Beem factor can be written as: 
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2%? W24525) (3.20) 
and 


My M, 


=— iN - 
mee, ) = UCzZ,,2,) + : » a aap M2 W(z,5Z,) 
m,=0 M5=0 . 


C321) 
(m,+m, 70) 
The inverse Z-transforms of eqs (3.20) and (3.21) are 
oem OF 
o Clee — B w(n,-% ,n,-2,) (3.22) 
Ve (lhe 17% My -h 
LO. Wet 
and 
My We 
w(ny,n,) = u(ni,n2) + iy » Sm wW(nj-M2,ni-M2 ) 
ule, 
Pr =O (2223) 
(m,+m, #0) 


Without loss of generality, Ls 1s chosen equal to M; for all 


1. The structure realizing eqs (3.22) and (3.23) is shown 


in Figure 3.9. It is noted that the number of delays in 


Meee 3.9 cquals the order of H(z,,z,). The structure is, 


therefore, by definition, canonic. 
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Merinition 3: 


The state Ang SMM 2? is defined to be a number which 
Sis 


Bee eoecethe time instant (non) the (pt1)5*, Catal) delay 


wa mie | 


3 direction, respectively. Following definition 


3 the states are marked in Figure 3.14. Comparing eqs (3.22) 
and (3.23) with the structure in Figure 3.16 easily identi- 


fies the following equalities: 
Xo q (myn) = w(n,,n,) 
X59 (Myon, ) = w(n,-1, n,) 


Xu 0 Grime) an ew Gite noe Tie) 


Xyq (nyon,) = CU p nigel 


X,,(n,.n,) = wnj-1,n,-1) 


M1 Gis) Te (GL IPRS gles cd 8). 


Xo, Gilani Gage tee vee 


XM, OR 5 Hel) = ic SLE eS 


“MM(ny3np) = w(ny-M,,n,—My) 
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From Figure 3.16 the relationship between the state vari- 


ables is identified as: 
Xoq (n,n, ) = O16 Xyg (mi-1 ne) + Ao X19 (i-1n2) eee 


+ Oa) Xp_g(m-l,n,) + O59, Xp Cn, on,-l) t+... 


30 


oe {0 


ll Xgz n,-1,n,) + O54 X,,(n,-1,n,) joe 


+ Gon BOG) eae + Qin Xow (m-lon,) +... 


+ om Xam, Cet abi eins 


X19 (ny, n,) = Xyqgfn-in,) 


e 


X97 (npn) = Xoo Cie = 1) 


X17 (n,n) = a ith 


Si sil) 


Kom SP aeMe) = Xow -1) 2 


Xiu (ny oN.) = ayy, Soe 
These equations are written in matrix form as follows: 


XGoeme ee Ae XCn longs Be X(njn,-1) + Bu (n,n) 


(3.24) 
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00 
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@ is °tme Zero Mmercrix. 


and 
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The output equation is similarly formulated in matrix form: 


ap eameme ccna np tec x(a j-1l) + Dutnon,) (3.25) 


where X is defined as before, and 
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2. State equation in N-dimensions 


The procedure for obtaining state equations in two 


dimensions can quite readily be extended to N-dimensions, 
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i.e., the equations for y(n) and w(n) are 


VCE = ) 85 w(n-2) (3.26) 
g 
(G<2 ie) 
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The state and output equations 
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ee 


D. DISCUSSION OF CANONIC SECTIONS IN N-DIMENSIONS 


ft was stated in Section B of this chapter that a one- 
dimensional high-order transfer function can be implemented 
in general as a parallel, cascade or hybrid arrangement of 
second-order sections. For reasons of noise generation, 
coefficient sensitivity, economics, logistics and others, 
these low-order sections are very important. The canonic 


Healizattons, for example in Fig 3.3, which constitute 
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the subclass of second-order systems having minimum number 
of components, are widely used today [2], [40], and [41]. 

The non-factorability of multi-variable polynomials and 
the fact that no stability conditions except for the two- 
dimensional bilinear case (eq (2.2)) were available up to now 
explains why recursive digital filter design has not been 
achieved in two- or higher-dimensions. 

In Chapter X, N-dimensional time domain design techniques 
are presented which result in high-order N-dimensional trans- 
fer functions or a parallel arrangement of low-order sections. 
The order of these parallel sections depends on stability 
and design considerations. It was shown that factorization 
Seepemtial fraction expansion cannot in general be achieved. 
However, approximations (See Chapter IX) are possible, so 
that the consideration of canonic sections is of value. In 
Appendix F, a state space approach is proposed to derive 
multi-dimensional canonic structures for the general second- 
order two-dimensional case. The derivations of canonical 
sections in general are not further pursued, but are called 
to the attention of future researchers in the multi-dimen- 


Sional digital filter area. 
Eee oUMMARY 


In this chapter an N-dimensional structure for the imple- 
mentation of digital filters was presented. The direct form 
Peewhich combines features of the direct forms 1 and 2, was 


used to develop compact realization of N-dimensional filter. 
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The significance of these N-dimensional realization schemes 
stems from the fact that realization schemes can be used 
directly to implement hardware, where components of the 
realization are replaced by summers, multipliers, and shift 


registers. 


Bae 
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Fig. 5.14: BLOCK DIAGRAM REPRESENTATION OF A 
SEPARABLE SYSTEM. 
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fie hi eae tnOD TOMDETERMINE STABILETY OF 
AN N-DIMENSIONAL RECURSIVE FILTER 


As INTRODUCTION 


The approximation part in the design of recursive filters 
1s carried out by choosing filter coefficients to approximate 
a given specification in time or frequency domain and evalu- 
mime the stability of the resulting filter transfer func- 
tion. The generalization of one-dimensional design tech- 
niques for recursive digital filters, in particular stability 
analysis, to the multivariable case 1s not possible because 
the fundamental theorem of algebra which permits one to in- 
vestigate stability of higher-order systems by applying 
Simple conditions to the roots of the factored form repre- 
sentation, cannot be extended to several variables. 

[Time domain design techniques are developed 
in Chapters VIII, IX, and X. At this point, however, it is 
of interest to note that the stability investigation 1s a 
decisive part of the design procedure. From the design 
engineer's point of view, it is very impractical if the 
stability of each individual filter design variation must be 
determined using complicated methods which result in time 
consuming computer use, and when the filter coefficients 
causing instability cannot be identified and, consequently, 
no knowledge is obtained to what extent coefficients must be 
changed to stabilize an unstable system. Due to this fact, 


the design of multi-dimensional recursive filter has not 
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been achieved. In fact, there is not even a simple design 
Method available [2]. 

In the following Chapters, an approach to the solution 
of the stability problem is presented which is not based on 
the N-dimensional roots of the characteristic equation. In- 
stead, it employs the N-dimensional Taylor series expansion 
of a transfer -function whose coefficients are values of the 


unit sample response, h(n), i.e., the Taylor expansion of 


H(z) = ) h(n) z (4.1) 


(<n. <=) 


H(z) is: 


It 1s known that an N-dimensional system is stable if every 
bounded input produces a bounded output (= BIBO stability). 
For linear, shift-invariant recursive digital filter a nec- 


essary and sufficient condition for BIBO stability is 


oes ) ini |e <a Cu) 


Proof: See Appendix C. 


To identify stability conditions or conditions for abso- 
lute summability of the N-variable sequence of eq (4.2) the 
general term h(n) must be identifiable. Four formulas are 
presented in this chapter for the computation of h(n) which 
are based on derivative, recursive, nonrecursive, and non- 


recursive combinatorial techniques. Both nonrecursive 
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formulas express h(n) explicitly in terms of the N-variable 
characteristic equation coefficients. Several examples and 
a comparative evaluation of the formulas demonstrate the 
perspicuity of the nonrecursive combinatorial method. 

The above computational techniques are employed in Chap- 
ter V to compute unit sample responses for the bilinear third 
and fourth-order (uncoupled) two-dimensional transfer func- 
Mfems. By application of derivative operator and Ultraspher- 
ical polynomial methods, necessary and sufficient as well as 
sufficient stability conditions are identified for the above 
cases. The concept of a "stability region in the coeffici- 
ent space" is introduced and pictures of the bilinear filter 
Stability region are shown. 

In Chapters VI and VII, sufficient conditions and neces- 
Sary conditions for the general N-dimensional filter, and 
methods to compute necessary and sufficient conditions for 
the general N-dimensional transfer function are derived. 
several examples are presented to demonstrate the application 


of the above methods. 


B. TAYLOR EXPANSION FORMULA IN N-VARIABLES 


The Taylor series expansion formula for a function, 


lage Of N-variables is given in [52] as: 


a . = 9 5 k 
ECZ) = >. fr (ak cee aN im) P(z)} Cue 3) 





69 











where 





denotes a differential operator. Since the expansion will 
always be taken about the origin, we could also term it 
Mc Laurin expansion. Following [53], eq (4.3) will now be 


reformulated into a more intelligible form. 


The powers of any term, 


Seales! , 2 _\N 
1 in a I Zn Choe 


of the expanded form of the product in eq (4.5) will satisfy 





“i Ly = k. The number of times the term of eq (4.6) 


occurs, equals the number cf different permutations of the 


% 
laeeaetor. iitlhus, 











V\ | ( k (: 5 ) (: 5 )" (4.7) 
. ; 1 5 eee N IZ 
ai ee a N 


mere i = [i,,..., ind assume all positive values satisfying 


Gbemrelation 1,... + Ly = k. 
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Peation (4.3) is rewritten, using eq (4.7), as:* 


_ ». » , al Oe 
4 = ; a a a a ar ar = 1 ° 
Key f gt rani i a BCzZ) Z 


Zz 4+. 821 () 


Gite. *. »ti,=k) 


which is equivalent to 


; »y ia ctamn : ~ 
F(z) = SEE TETE ‘ : rae TAS | Zz (4.9) 
eee oN ues aN 





Also, if we define 





: VE a2 aes er 2 
fo ai ea) fen) elle C7) (4.10) 
9z,7 Z., N 
oN 

then eq (4.9) becomes 

ez) = » pS pote bee eu joe (CE, ak) 

1 eee N° Ci) 
a 

CO Snes) 


This form of the Taylor expansion of a function F(z) will be 
used in the subsequent derivations in this chapter. Its 
application will be demonstrated in the following short 


example. 


Example SA lies 


The first few terms of the Taylor expansion of a third- 


order two-dimensional transfer function, 


¥% ° al ba N 
N lass tie, eh N ? = = 3 = 
ote Y) (ars) (2y ye) FC) = y) z a sai iy ) F(Z) 
¢ Zz 
N 








2 Chia) 
Lo ~ G19 Zy- %, 2 -%q 2 


can be expressed, using eq (4.11), as: 


1 _ . oe 
BEC ZisZ >) = 3 yom Lt boty 2425) | Zee 
1 SG 


: (0,0) 
1, =0 














aa) 7 
jy 2 eae ) Z 
1 ; 0,0) OZ, ,0) 
> 8 Bos 
+ GZ Zo) Vb, Te ~—— —, L(z,,z,) | zz, + 
OZ, 92> 2 an 172 OZ) aA 1? "2 0,0) 2 
= 1 + 19 2, F Ora at eng 2 2 + 
2 2 
+ (3a, ; By PS Chon Say) tap AOS (4,12a) 


© GENERAL TERM OF THE N-DIMENSIONAL TAYLOR EXPANSION 


It is useful for later derivations to formulate the Tay- 


lor expansion of the N-variable transfer function H(z), i.e., 


H(z) = ————_— (4.13) 
io OG? 


moeterms of Zz This will be possible if the expansion of 


ae i 
L(z) = To te (4.15) 


can be found in terms of powers of 2z 


Te 





iaeearee Of the following sections it will be shown that the 
coefficients of the Taylor expansion of L(z) can be explic- 


fee derived in terms of the Taylor coefficients of Q(z). 


ime Derivarive Formula 


If the Taylor expansion formula of H(z) is written 


as 


Cz ee Y RG yet (4.16) 





h(i) = = z all gy (dies + in)zy | Claes) 
ies * N° (0) 
2. Recursive Formula 
Theorem 4.1: If the Taylor expansion of Q has the form 
ay = a,z. zt : ee Cie) 
} (i) (5) 
a 
(O<1i,<~) 


and the Taylor expansion of 
~“< 7 il 
L(z) = Cisal ch) 


ia CZ) 


ics 


ie 





Sar Y) = al 
then 


en) = ) G5), OGD ey 


Proof: It follows directly from eq (4.15) that 


CA) SG). Claes wale 


eipestituting eqs (4.18) and (4.20), it follows that 


7 1 = i r | i 
7 


b+ 


ae 


When i = 0, since 25) = 0 it follows Did) sek, cree eall 


feeiies of 1 £ 0 it follows that 
io. Oo I if 
>») va) a= >) bea) Z »» any 


As a result 


°(Z) © ». “i-9) Pj) 


j<i 
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Fxample Ee 


The first few coefficients of the Taylor expansion for 
L(z,5Z,)> where 


il 


eee es es Cue?) 
Lo Oy 4~ Oy Ze - GY 


Mi 7A, A 


are computed by eq (4.21), where the only nonzero coeffici- 


ents of Q(z,,z.), which is defined by eq (4.18), are: 


ems 1° © F0,1) ~ B98 Bean) ~ OP? 
al 1) 
et i ) = a y 2a . ° b e e 
3 2 oe ae Gaia 2-0 2.) (J Lodi ad 
Wo pe 


jJitje *¥ iitie 


itt ol ilOws: 


ile an 
b = : a l=a 
Co) \ ). Coa laab 10 


a i a 
J, vie 


1 
b _ = a : ° = 
(0,1) d. ) tet h? 2 53 
420 7,500 


ih Mig 


15 





Di s 
YS » oe Ge =e baa) = 20, , ig 


tied ane 
eat ab #1 


2 


cae ) 


al ; 
CZs i) = 85g? Oy, 


iar 


2 1 
-——4 
(2,1) : » y Ses o> Gyeai a 


> 8h 2G) | Sano) 8 ean 


These coefficients are, of course, the same as in eq (4.12a). 


The computations for the Dee 1) coefficient can be dis- 
>) 
Played graphically in the following manner: if the coeffici- 
S@eemor the 2-variable characteristic equation of L(z,,z,) 


are arranged in matrix form, as follows, 
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and if we arrange the b y in matrix form, 


Sdinos Ty 





then the equation defining re, ik 
3 


: 1 
D(2.1) i » = “(2-jp 1-4) hdd 
=0 


=) 


2 
3 


WM CU. 


can be visualized as the rotation of the first matrix by 

180° and the superimposition on the second matrix, such that 

Yo = 1 1s above Di, iy: By summing the products of stacked 
3 


coefficients, as shown in Figure 4, the result is: 


De) = S50 Dens 4 O61 (2 0) ‘ G29 7 (0, 1) 
+ + 
et a) 


eee Nonrecursive Formula 


To compute the term b;z, using the recursive formula 
requires the computation of all previous terms, i-e., Ge 
where } <i. To circumvent this disadvantage, a nonrecursive 


formula is derived in the following paragraphs in which 


Biz) 1s uniquely computed using the set of coefficients in 


OZ) . 


iif 





Q2 3 


O13 


O03 


Q22 


Q12 


9 2 





fier: GRAPHICAL INTERPRETATION OF EQ (4.21) FOR THE 


be 4) 


COMPUTATION OF EXAMPLE 4.2 
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Equation (4.15) can be rewritten as 


cz ey , | cz) ]" » Q*(z) Cua Ze) 


R= 0 L=0 


Then, if we define 


M,(Z) = Q'(2) (4. 214) 


eq (4.23) becomes 


Le @) 


ez) = yy M, (z) (4.25) 


L=0 


P@emecicine the (1,,...,1,) derivative with respect to 


N 


(Z,5---52,) of both sides in eq (4.25) leads to 


N 


potis-++2tN) (3) 7) mM, tne 22N) (2) (4.26) 


&=0 


where the notation of eq (4.10) is used. 
Mt is concluded that the Taylor coefficients for L(z) 
expansion can be computed whenever the term M Suipe scone) 


can be calculated. 


Theorem 4.2: 


ie (x ) = n(x)...r (x), where r.(x)-isa function as 


one variable, then as shown in Chrystal [54] 


19 





Theorem 4.2 applied to eq (4.24) results in: 


(0. 55e0 oe A, Go Oe 
oe. OT a ie fae . =e 


ati oe cc 5) (ly. 28) 


if Theorem 4.2 is applied again to eq (4.28): 


iM, CLartar0o-+ +103) >. ». 


Jiohe ++ Jig=21  Jact---Jag=t2 


9 T102J32029 +++ 9I(Z) mg Fige Jog 295+ - 9X2) 


(4.29) 
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Then by induction, it follows directly: 


M, Stir++ + tng) = ~ oan yy 


oy eae Inpt +++ +Ing =1y 


rT ly 
Joss ++ stig InNooc ++ Ing 


Fin - — >» IND (5) ove gSIip> ae > INL) (Zz) 


Ca 309) 


Combining the coefficient term of eq (4.9) and eqs (4.26) 
and (4.30): 


1 Cae ee 
.3,? * N“(@) -) » 
a=0 [3 


hats uly vee 





) 


I not sets INg=iy 


il ure Ll 
Jiot +++ IQ! Ino Iyg! 


gSIio- +> 2Iyo?... g(In>-- ° sINQ? 





(4.31) 
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EMoyanes (4.16), (4.17) to 
O(z) = ) az) ye Se Pele) = 0 
us 1, <=) 
meeecocificient term (7) is computed as 


1 (iy,---,4n 5, | (4.32) 


eel) ales! ° a) 





Similarly, for eq (4.20), 
Ee > ee 
mic 7) = Dia) a (4.33) 


where 


(1 peer neal: ) — 
Bi) eit 4 oat L N CZ) | (4.34) 
1 j ea die: N° (0) 





Paembay eq (4.31), DOG) equals 


il al 
doc “Jig! IN’ -+Iyg' 


o Fro? >In? ae gat > ue? (0) ee 


ieeeg (4.32) is applied to eq (4.35), then the result is 
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Summarized as follows: 


Theorem 4.3: 


If the Taylor expansion of Q(z) is 


a y ee: 


(O<i, <@) 


(436) 


and if the Taylor expansion of L 


7 | 3 
L(z) = Y) biZy 2 (4.37) 


Jil mee dy yecn J no 


Ge), 


Bil aoe ; 
Chios ++ > INQ Chap os + + Inyo? 


Using vector index notation, eq (4.38) becomes 


Cae » Y aa eee (ing) (4.39) 


where the multiple summation is performed over all %#-tuples 
of vector indices Gyo 29 again which sum to the vector index 
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[meine application of Theorem 4.3 is illustrated in the 


following examples: 


Example ut aie 


For the two-dimensional transfer function 


: | 
L(z,z ) = —=————_+—__—_ (4.40) 
w ie Sh Sty Se Coie 


which was also used in the previous examples,the first two 
coefficients of the corresponding Taylor expansion are com- 


puted as follows: 


L=0 Jiot- ++ try Jact- ++ *]a,g 


ary-s : cee Are . 
(ak co dane C31 Joy? 


where the only nonzero coefficients are 


=Gao)  » 770 (4.40a) 
Big1) = %, C4.40b) 
Foon) = M20 (4.40c) 
Mores = © the indices j:p and j25 Sum to i, = 1 and i, = 0 
Wiemever j,, = 1 and j,, = 0. 
For 2 = 1, the indices sum to ee ed and i, = 0 whenever 


CJiodpp and (429.32) belong to the following set of ordered 
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pairs: 


Grohy = 40.195 G,0b and G,yi,) = {co,o} 
Hee. , 


1,0 Lety {5,¢255, 70}, and 5,905, 79) 


For & = 2, following the same line of argument, 


Sies dq Jy 2 {(0,0,1)5(0,1,0)3(1,0,0)$ and 


Ggrdy edgy) = 4€050,0)f 


fever the indices for 2 = 0,1,2,..., and the coefficients 


of eq (4.40) are substituted in 


aa) y y SS pda) (2 = 0) 
fl <9 

+ ee ms R= 1 

\" yy mG ee, Giri) 


Duele Jo oto. 


=1 =0 


+ . , : . ary: : 
a » “Gh coda “(Gh se dae (F129 J22) 


Jil thi. Joot Joi Jo2 


=i =1(9 (2 = 2) 


then the equation on the following page is obtained: 
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Se 


Cen : ee 0) 


2 0) Gay BO NeG.e) 4 
+ A a a rae 2 a 

J) IS OE Stan» COR OL So) Ga) 3s iy 
+ 


Gays Con Coo) 


= Q19 , since a = 0% 


(one 


Similarly, 


Se) = 2604 


meeomple 4.4: 


The term, Dos 1) 1s computed using Theorem 4.3 for the 
3 
feterer Tunction i(z,,2,) , which is defined in examples 4.3 


as follows: 


Hemme = 2, 1, = 5, and 1, =1, equation (4.38) becomes: 
b = ne 
eset) > . . Niele Gi s'J2 ) 
2=0 Jig Wghse 
25 =1 (4.41) 
or 
+ (2 = 0) 


en (58) 


86 





ep) 


Inge dhe 


ans Ae 3 Sey Coe 
iy > lanpdan Ahad Sina 2 


Joot Jai 
+, = = 


Jioth 1 


+3,25 


Clueas ‘ 
Ad. » only 


Jog 1, 1 da 
+4, sae 


ao LJ wy 
t= = 


eee 
(3139 b 3 


o?. 


Jiro Jui Jy2 


ase ° 
yy C310? Joo? 


Sa See 
43 0 Joy Jo2 


1 he eae 


ars 
C14? J, 


BR OD 


oan 2 


Aes ; 
rr J, 5) 


a,- : 
Sane ee 
Jo me WJ, 2 


Jnecanonele ets 


=5 =1 
Aes : alicnas ; 
eee Sarasa 
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See «ie 
Nioston. 


“Ghetoe 


SG = 4) 


Eee 
Guns ‘ Arye. F 
Sapo da? Shae lage 

C= Gen 
“Cio y Bch 

Cae 
aio : aps : 
Seda ree) 

(C= 15.) 





The only nonzero coefficients in the characteristic equation 


Pik 2 ,2,) are: 


(Gl) = As cl) 
= (0), |) an amor | (4.43) 
Serio) = a5 x20 (4. Un) 


It follows from eq (4.41) that the sets of indices 


{a1o -++sd:,} and eee, 


must add up to 5 and 1, and the set of ordered pairs: 
{aed Os Gigrdy)} 


must correspond to nonzero coefficients, 1.e., as listed in 
eqs (4.42) - (4.44). 
Applying this reasoning to term (%= 1) the two sets of 


mieme@es are obtained. For 


Ghote) > {(0,5)5¢1,4)5(2,3)5€3,2)5(441)5(5,0)f 


ema for 


Gres...) ; 4(0,1)3(1,09$ 


It can easily be verified that the products of coefficients, 


which is defined by term (2= 1), becomes 


88 


> > fal Geeta?) Mame Gora 


J, 0°35, ‘a Jo pt Jo = 


r *( 0,0) COR So) 2A sai) F(2,0) a Gsest) 


MECN G.t) = «4 (4,0) 7G,1) = ,0) 2.) 


MGs). Gt) 2G. ey AG) Gee) 4 


350) “OR em) Sein) 2 aren) Sid =e 


Following the same reasoning, the term (2= 2) becomes 


a;-« : a;- : a : _ 
d. » ©3102 J20? CHip hp (J, 25,2 = 8 


jis dae pin 5% 
=5 =] 
The two sets of indices for the term (2 = 3), which add 


to 59 and 1, as well as corresponding to nonzero coefficients 


in eqs (4.42) - (4.44), lead to 


ae Eeareaeee era. Alaa 
q >. “Sapp den  SInnddawe ihes ae Cues ie, 


Jno Toit hat dee 
is Joa = 


= Meso) (2,0) ~G,0) “(041 ) ¥ (On): Ss aon “Bay” 


STIG Ce Gn) 8G) 2,0) G50) 85a) 
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emmeino) (no) (.1)s (1,0) a@no) -(, 1) 72,5) 


(20) Kara) ea Ca 0) 7G 0) 4 0,1) Sea) ace 0) Co) 


+ + 
Mao) Gn) Gn) Ca) COMO =OMo Corie cs GEO) 
cremac ae) *(ime) “(mmm = F02,0) 20,0 va ,1)7GIo) - 
= l2 ce O19 Ay 
Similarly, the term (2= 4) is 
J >. = Giesao ns Gunes) (Ghpipe) Grenieo) Ginestae) 


Ji of I, the Joot Ie IIo 2 


Jha 9 J #I,,21 


= 3 
cee G15 #80, o,, 


and the term (2= 5) becomes 


a;-. ee. ee : ally : a : 
4 ) Saige lee STpee ds se Te daw Sipdad 


ot Hh 2 ie Fee 


14 ee Ine dao 2! 


a - ° Gy. 5 
Bry 4) See lee 


It can be shown that the terms for 2 > 5 are identically 


Zero. 
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nus : 


2 2 3 5 
Dest) = 12 a), 4, 4,5 + 20 Gy yy Aq + 6 Og Hy 


these examples show that the application of Theorem 4.3 
is straightforward, but not trivial. The computational 
method allows the determination of some b without requiring 
the knowledge of any previous b's. Tee Ooms section 
develops a nonrecursive combinatorial formula which has 
several advantages as compared with the above outlined tech- 


nique. 
u. Nonrecursive Combinatorial Formula 


Following [53] and using the vector notation defined 


in Appendix A, l.e., 


SP Glens) 


it is defined that two £-tuples 
Cyoveeodg) and Ciys-- 5) 
are equivalent, i.e., 
Slope 0a9dae me (Giigescoda } 
whenever any given vector index j appears the same number of 


times in each &-tuple. Therefore, if 


(Gype ody) 9a Gl ape saeay ie 
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then it can be easily shown that 


= CE aia ay ) = ay 


i eG) 


Also given an &£-tuple Gece of vector indices for each 


vector index A tS defines the number of times the vector 


index 5 repeats among (Flea ooo Gay)! ier follows Gdirecully that 
os 

a . ececsvo a 7 a a aoe C4Y.45) 
er *) ee eb, 


For example, consider the term 


= (iy) maaan Gap = 20) @ Ge i) “Or Cp = Ole. 
where = = {2,95 (0,295 (2,0)$ 
3 
on co Pty ENG) Cts al OU al One (Can 


It is also observed that the sum of tS is equal to 2 + 1, 


» os See iad | Chita) 


d 


se. , 


If there are (Spy ++ 5S,) distinct vector indices among the 

&-tuple, Gipae sess): then et 4 (6) meio tetas erey eee ce Since 
Zz n oe 

another &£-tuple (Jyo- +09 Q) is equivalent to(j,,-.,j,) if 


and only if each index (Sio+ ++ 98,) appears alike 990 dea 9) 


the same number of times Us »-+-+,t5 it appears in (Tecacniy? 
0 


92 


it follows from combinatorial theory that the number of 
£-tuples equivalent to Glance 4 oo ig) Geweuls the multinomial co- 


efficient, 


t- + +t- \! 
Qt i 7 sats of BS ( 0 he 
e State EH | 
ts. ion “a ees 0 k 
CET) 
Note: 
Sy pe ewes etinedywilth k = 0,1,... for each 
distinct index (i peeesdy) which appears in the coefficients 
mG@(z). ce is a positive integer. 


Combining these results, the following theorem is obtained: 


Theorem 4.4: 


If the Taylor expansion of Q(z) is 


Z yy i : 
eGz) = a5) ZL 5, 20) - 0 


a 


(0<i,<@~) 
and the Taylor expansion of L = 1/1-Q is 


eae: \\ : if 
Gay = DOG) Z 


1 


(0<1, <~) 
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Example eye 


For the two-dimensional transfer function 
previous examples, the first two coefficients 


sponding Taylor series expansion are computed 


y |. 


Ce 0) (0, 1) 
awe 


 @. che 


a 
(1,095 (0,1)9" (2,0) C1, 


EO o SE Co 
OG, 4) 
10 


where tao) (1,0) feet CO od tz 0) (250) 


re) 


For a wo and a, = 0 


Gl 0) Fo, 1) 6021) + Ong) oo? 


"Ca,0) 


Therefore, ae = Gem) Oso 


Similarly, b <= 


(0,1) 01 
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C20) = 


CIE eke) 


UiSevel Glin) ele 
of the corre- 


as follows: 


0) “(o, 1) §(2, 0) 


) Shey Ie. OD 


ele) 


Cr) 





Example Wyo 


For the transfer function , L(z), as in example 4.3, 


compute b’, ,)-+ The only nonzero coefficients in Q(z) are 


=o) 10 
OG 
B(o.0) ~ “a0 
Thus, Sa aeGL A 0n 
Seean Ore) 
St (250) 
and Tao) 61,0) a Te20) 62,9) < to 1) 6051) ="@ 1) irom 


eq (4.48b), whenever 


ey = 2 eye eae Soa) 7 Me oe 
ern SGenye Snes oe OS 
Seay ee Eee = e tio) = 1. 


Therefore, substituting into (4.48a), 


Sel Slum leaatt 
Be. 1) = Sn Clone Cre ois Oyo Chg e 
el 2 Joi. 
Sy db 
u O19 My, = 
Sige! 


2 3 5 
12 yg My %yqQ + 20 Ayo M1Geo + 6 A, OH, 


2 





She Comparison Opeeeowlus 


The four computational schemes presented in Part C 
are summarized as follows: 

The derivative method is very difficult to apply 
Since the computation of the Dee coefficient of the Taylor 
expansion of L(z) requires taking the (ij,..-4,)) partial 
derivatives of L(z) and evaluate the result at zero. This 
technique cannot be easily implemented on a digital computer. 

The recursive method requires for the computation of 
ye) the knowledge of all Day where } <i. This technique 
is best suited for computer implementation, especially when 
Semputing a finite set of coefficients, i.e., {bz} , where 
it. sl.. It will be used extensively in Chapter V. 

The nonrecursive and nonrecursive combinatorial 
formulas relate the coefficient DZ) explicitly to the co- 


efficients of the N-dimensional characteristic equation. 


To compute Dex 0) for the transfer function, 
2 


le 
L(Z,,2,) - 1-oa..Zz - a ne 
10 “1 Qo, 42 20 “1 


Using the nonrecursive technique of Theorem 4.3 requires: 


(k-mtj)! 
“>! m! (k-2m)! 
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multiplications in addition to identifying the same number 
@amindex sets, where [X] means the largest integer < X. 


When applying the nonrecursive combinatorial method 


of Theorem 4.4, only E + 1 multiplications must be per- 


formed in addition to identifying | + | number of t's and 


mieetinomial coefficients, which are tabulated in [55]. The 
following table further elucidates the comparison of the 


number of multiplications to obtain Dy ) for the above 


kJ 
mes Zz (+): 

@4)°° term Number of computations 

k = } = Th 4.3 Th 4.4 
0 0 1 3h 
i ab 2 if 
2 Z S] 2 
3 3 a) 2 
Uy + 0 3 
5 5 924 3 

ig) 10 3045185 6 

LS) 25 ese See ya aU 8 

20 20 Uo 34 Eb + 13 ae 
25 25 1.881 E+17 | 13 
510 30 OEE ci SI 12 7a ae ie pee 416) is 


The tremendous savings in computations when apply- 
ing the nonrecursive combinatorial method of Theorem 4.4 are 
obvious. 

There 1S one important application of Theorem 4.4 


which will be used in Chapter VII to establish a mathematical 


out 





formulation for the propagation of transfer function coef- 
Mietents in unit sample responses in any dimension and any 
order. These rules are then used to develop in Chaptens vii t 
a transfer function extraction algorithm, which in turn is 
used in Chapter IX to develop time domain design 

techniques. Therefore , the high Speed characteristic 
which evolves from the few multiplications is exceptional 
but the mathematical form is of greater importance for the 


Subsequent chapters of this thesis. 
PUNT T SAMPLE RESPONSE 


By comparison of eqs (4.13) and Cro) sat he relat 1Omeie. 
hez) = ez) is eccculvetdentitied. The following theorem 
establishes the relationship of the unit Sample response 


entry h(i) and b _ ,; where h(i) is also the 7th Geer ficient 
(a) 
of the Taylor expression of H(Z). 


m@eorem 4.5: 
If the Taylor expansion of OC2) wis 
=. : i _ 
OCz )! = ariy 2 : 2 (9) 0 C449) 
a 
(Osi; <@) 


and if the Taylor expansion of ez) = I SGy is 


Lia) = » by) zt (4.50) 
iL 





J. 


then for the Taylor expansion of H(z) = ee in the 
form ye Ce 
H(z) = ) GGL) (4.51) 
1 
(0<i,<«) 
each and every Deas = Bee) ‘toe esl dLaL i, and the coefficients 


n(i) are related to a(j) as stated in eq (4.17) and in 


Theorems 4.1, 4.3, and 4.4. 


Note that the notation h(i), instead of h(Zy> 1s the 


ea 


commonly accepted one to designate the i coefficient of 


the unit sample response. 


Table 4 summarizes the h(i) versus 2 (7) relationships. 
Equation (4.51) can be written for the two-dimensional 


case aS matrix product: 


step) Ze ze oe) A Cl 2,7 257s.) 


where 
h(0,0) h(O0,1) h(0,2) h(0,3) 
COE KCl) nine? ) lena), 3) 
7 OD nn 1 G2 nn VG > mE VG 
i 
= NEON ieee) hGaee) nAcssa> 


ome n(0,0) = 1. 
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Be Cg Sapo eee ai) 
“Ie c @S 
U ; Al eee = 2h 
> 0 Oe kee 
S 4 24 
(2s eee jen Ae 0c = (u)Yy 
24 ey Bate + 29 
Beer ORC: % 
% 
Ce. OL ( va (u)U 
aN 
u>e 
(Ey CESU)e « = (u)y 
(0) Nao... qtu 
Rieger: T eo 
( oe Use) 
eTNwto0g 


Cir tr) 


(30) 


(T° i) 


(Wat09U4 ) 


Shh 


BE fi 


ee an 


eae 


uotqenby 


(Z)H NOILONO4S USASNVHYL NSAID V WOHd (U)U SHTYLNG 
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V. APPLICATION OF SERIES METHOD FOR THE DERIVATION 
OF STABILITY CONDITIONS IN _TWO-DIMENSIONS 
SL  PEMENS TONS 


[a =LNTRODUCTION 

ine formulation developed in Chapter IV is applied to 
compute the unit sample response matrix {h(€n,,n,)} of the 
third and fourth order uncoupled, as well as second order 
coupled transfer functions. Conditions for absolute 


mee ecnce of matrix {h(n,,n,)}, 1.e., 


2 = yg > [n(ny,n,)] < © 5 21) 


n,=U n,=0 


which are proven in Appendix C to be identical to BIE© 
Seability conditions, are found in two steps: 

(a) for the partial infinite sum of absolute entries 
in each column (row), and 

(b) for the infinite sum of column (row) Sums, where 
the sequence of Summation steps corresponds to Cauchy's 
method one (two) respectively. It is shown in Zon Gl esc ees 
that the summation of Successive diagonals, which is Cauchy “s 
method three, leads to the Same result. 

For the third and fourth order uncoupled transfen FUNG Eom 
the latter method is used in conjunction with Daeine@mare I theory, 
resulting in a sequence of geometric series for which 
convergency conditions are known. In Cede tOeap ply. the 
binomial theorem in forming the sum of successive diagonals 


in {]h(ni,n2) |}, the absolute signs must be ISSHNC etal WW i@alte\ ais 


iG) 





ieiepeexcept tor the trivial case, by identifying each entry 
Min; .n,) to represent an Ultraspherical polynomial. By 
using its known characteristics, necessary and sufficient, 
Seewell as sufficient stability conditions, are derived 

(see tables 5.la, 5.1b). 

It 1s also shown that in the second order coupled case, 
Cauchy's methods ae or two in conjunction with derivative 
operator theory results in a sequence of geometric series 
for which necessary conditions for column convergence and 
necessary and sufficient conditions for column sum convergence 


eee adcentified (see table 5.1c). 


Pee COUPLED FILTER 
Definition:  Uncoupled Filter 
A two-dimensional digital transfer function is uncoupled, 
if its two-dimensional characteristic equation has no cross 
Say 


terms, 1.e., terms of the form omy Ze ae 
] 


Tg 
No : 


where n, and 


No are both non zero. 


The general uncoupled two-dimensional transfer function 


forme De written as 


: iu 
CZ) 5 Z2) ~ £(zi) - hizo) Coe 
= y [h(z2)]P teleay}oe (5.33) 
p=0 
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YHLTIIA (AH1dNOO) HACHO ANOOES You 


Misia) TE ae S- | Tt0-T 00 | 


0ln —- T > | TT o4T O10 | 


>: ATESSSOON 





lover otstoya Malm cre ays, 


SNOLISNYWNIG-OML NI 
SNOILIGNOD ALITIGVLS ‘T'S Wave 


zz t7ttlHy - C7t0H = 


(@¢z°'z)y uofFqgoung dejsuedy 


IGS) 


iMmelnara Order Uneoupled Transfer Function 


For g(zi) =1- Are ig - ee Cony 


h(€z2) Sin ee C5235) 


Peco. 2) and (5.3) become 


L 


(5.6) 
= 7 z 
ji = AeA = 020Z1 - re 


CZ 15Z2) 


~? =«(pt} = 
1 ) (p ab yzor 


= Il 
(1 - Q®1921 = "8902 


a0 (Se) 
Equation (5.6) is realized in direct form as shown in Fig 
5.1. Note that the direct forms are identical whenever the 
numerator polynomial is equal to unity. Equation (5.7) can 
be used to implement the third order uncoupled recursive 
filter in hybrid form, i.e., non-recursive structure with 
recursively realized components as shown in Fig 5.2. 
a. Unit Sample Response 

The planar unit sample response A of the two- 
dimensional third order transfer function (5.6) is listed 
in Table 5.2. The details of computations, which are 
performed using the recursive formula of Theorem 4.1 are 


contained in Appendix I, Section A. 
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a 
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Fig. 5.1: DIRECT FORM REALIZATION OF THIRD ORDER 
UNCOUPLED FILTER IN TWO-DIMENSIONS. 
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y(n, in.) 


X (n,n) 
. pt - 27-0. 25 1p) ——h 
i0 1 201 : 
-| 72 
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a. © 
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- | oo 
Zi 2] 
a Vy 
O01 -2-3 
(l-a Z-a Z) A> 
0 20 1 A 









1 


- 2] 
M(l-a@a Z-a Z) 
10 1 #201 





ma. o.c* HYBRID STRUCTURE OF THIRD- 
ORDER UNCOUPLED FILTER IN 


TWO-DIMENSIONS. 
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be Summation of PpcOmmpe Entries of 
Unit Sample Response A 


The sum of all absolute entries of column 0 








becomes: 
oe) = J |h(n, ,0) | 
n,=0 
oe [n1/2] = 
(njy-m)! of 2m Sesh Cs 8) 
m! (n)-2m) ! : 
n,=0 m=0 
Semilarly, 
a [n,/2] | 
-~2m m 
(ny-mtl)! gi} 
R- ) | (5.9) 
n,=0 m=0 
and for the oe column 
00 fae | ! 
= : -2m j m 
a Given aon O, Q 
5 , » » 7! m! CRED E | Soo dll, 
n,=0 Tite) 


The sum of all absolute entries of the unit pulse response 


S = 3 a Jn€n,,n,) | 


n2=0 n,=0 


1s now easily identified as: 


iDia — A 1G m 
+ - 1 1 2 
mm) ) |) Sere can 
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ec. Necessary and Sufficient Stability Conditions 
G10 
m 





OIG CG using Ultraspherical 


Polynomial Method 

It is proven in this section that the two- 
dimensional third order transfer function (5.6) is BIBO 
stable if and only if 


Qeo < l- loro | - loo i | Coy LD: 


Q 
mm RN 
fo) 





Pimemever Goo > - byetnspeecttoneor eq (5.11) efinst for 


pI ay 


miemcase O29 > 0 and secondly by utilizing the known 
characteristics of Ultraspherical polynomials. 

GP etchant, Meochadttrons ror G55 2 0... It as 
proven in Appendix C that a two-dimensional system is BIBO 


stable = aig only if 


> ) lh(€ni,n2)| < @ C 5:ie3) 


i SC igi 16. 


Equation (5.11) can be rewritten as: 








fia 7 2 | 
a va 1 ie (nyt+n2-m)! O12 9 ks (5.14) 
h2=0 n,=0 m=0 


Laue 





which for Q25) > 0, becomes 


[n,/?] m 
Ne (nitn2-m)! 29 
4 3 Jara | oo’ Tg! m! (ny-2m)! ( S24 \etinad 


Mg=0 ny ,=0 =0 


Since O19 and Qo; are real numbers. 
The summation of all entries of <A> , 
where < > is defined as taking the absolute value of each 


entry, 1.e., the Mea) term; using eq (5.15a) is 


[n,/2 m 
n Gor ramen) 
h(n; sn2) | os loro] Nex] z Soe (33 | (5.156) 
m=0 


1s performed following Cauchy's third method of summing all 
entries of an infinite by infinite matrix. (See Appendix K). 
With the detailed computation contained in Appendix I, it is 
noted that the sum of successive diagonals of < A >is 
greatly simplified by application of the binomial theorem 
and derivative operator methods. 

The summation of all entries of A 
converges or equivalently the third order filter transfer 


function H(z,,22) is stable if and only if 


O20 < 1 --laio| - [oor] Oro mCignes Os eeGoe ic 


See 





Notice that the validity of eq (5.16) is easily checked by 
inspecting the unit sample response matrix A as shown in 
Table 5.2. 

C2)> Stabalimevrcondrusens for 0 > aoe 2 = 
It is shown in this paragraph that eq (5.14) can be rewritten 
for a2) < 0 such that the summation inside the absolute signs 
can be identified as an Ultraspherical polynomial. In order 
to simplify the summation of all absolute entries of A by 
application of binomial and geometric series theory it is 
necessary to remove the absolute signs, which can be done 
utilizing the known characteristics of the Ultraspherical 
polynomials. 


hon eomanOemegm@s ld)ecan be rewritten (5.17) 


in the foliowing form: 











: 7 [ni/2] n,-2m 
$ = (Wage)? of a oe 
20 01 m! no! (n)-2m! 
= ¥—-20 
n2=0 n,=0 m=U (5.18) 


But the summation inside the absolute sign is identified 
using the theory presented in Appendix D, as Ultraspherical 


Memynomial, i.e., 





[n,/2] . n,-2m 
eet) mend. = : (ayy) (nitn2-m)! Q10 
my i » m! no! (nj-2m)! Jing 
Cor: 9)) 
m=0 


aes 





It is proven in Appendix D that the Ultraspherical polynomial 


Chie + 1) 0 
a) Se DOS tae 
ee 2¥-A20 
whenever 
Sk (5.20) 
Zi 280 
which implies a;, > 0 . Therefore, eq (5.18) can be 


rewritten as 


4 3 ny he Waal ae LL, | 
( =are | ote eG ie) Comers) 


1g 
n>=0 n,=0 : 2/-G2 0 





Also, for 


3 reel CS222)) 
2V¥-a20 


it is proven in Appendix D that 


C —__ O19 
ae 2V¥-020 


is positive/negative whenever n, is even/odd, respectively. 





Thus in this case, the absolute value of the Ultraspherical 


polynomial becomes equal to 
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2 (si)? a fad Ge Gy, 
Nn} 2V¥-A29 








Cn tl) 
C 419 
ey (- 7 


It is also proven in Appendix D that the right hand side of 


Bam 5.23) 1s equal to 


Gnio 1) a 
=e ——_+ 2 C5 7 8N) 
2V¥-a20 
fonda tion 
ee Ee ae (5.25) 
2¥-Q20 
implies a,, < 0. Thus 


Cnotl) (ntl) 
cl ees =C. —eswll 
yl 2¥-029 mt 2¥-A29 


and 


y y (ntl) 
CAG lees La (5.26) 


Oa ae 
No- =0 n3i= - O20 


eS) 





memautaons (5.17), (5.20) and (5.25) are now combined as 
0 > G29 > - =e (5.27) 


miomsimatiarly, eqs (5.21) and (5.26) to 


‘Gavia il)) 
= emotes” (lal 
-420 


ne=0 n,=0 


which by eq (5.19) 1s equal to 


[n1/2] 
m 
lox ie 9 {M2 (njt+n2-m)! Qo 9 
“ nn Gem! | eae 
n»=0 n4=0 — Cae Zee 
which is identical to eq (5.15). Therefore, it is concluded 


that H(z,,2Z2) as defined in eq (5.6) 1s stable if and only if 


ad2o0 < 1 - Jaro] - |aoil (5.29) 


Or 


(hic BO, 
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2 
(eect tlei1ent: Stability Conditions for a259 < = pelo 


It 1s proven in this section that the third 


order uncoupled transfer function is stable if 


lao1 | < (1 - v¥-a29 ‘a ox Qeqo < - Setiios Ce .31) 


by using the known characteristics of Ultraspherical 
Betynomial. For 


Q10 
S20e > see 


which is the condition for complex roots of the one- 


Mimensional second order transfer function: 


if 


1 - onan a ras Oi 5 


HCZ1) = 


it is shown in Appendix D that the Ultraspherical polynomial 


+ 
—" 1) O10 
ee 2V¥-Q20 


has ni real and distinct roots. This implies that the sum 
of all absolute entries of A must be performed without 
removal of the absolute signs enclosing the Ultraspherical 


polynomial part in: 


Aces 





Ny Nh 
5 = (Y-a20) |aoi| 


n,=0 n,=0 





5th) 
C_ e—_ (5.32) 
1 2¥-A9 6 / 


If we apply Cauchy's first method of summing all entries of 
an infinite by infinite matrix (See Appendix K) then the 


eum or all entries in column OG becomes 
CO 


Ny 
oy (Y-a2 0) 


n,=0 


oe Bona (5.33) 
; 2¥-a2 9 





It 1s obvious from Fig D.1 of Appendix D that there exists 





no range of (ts — ) where all om are of same sign. 
1 


2Y¥-a20 
Therefore, an approximation from Ref [55] is used, l.e., 








en tb) ee 2 esta 
Ch ee) < CSmst ) 
; 2¥-a29 ny 
Thus 
fo 0) 
ny iger nd. 
S4< (/-G00) oe 
ial - n, 
n= 
= pee ee. a 51 ee (5.35) 
a V-a20 Dy =) 
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and 


j 
me) s. eo ) (ie (5.36) 


which converges, if and only if 


2 


loor] < (1 - Vaz ? qvetdt 
e. Improved Sufficient Stability Conditions 
2 
Reve A296 < - = 


It is proven in Appendix I that a two-dimensional 


third order (uncoupled) transfer function is stable if 


@| 33 


sin 0*v-a29l 02 + In? ao] [1 - Kan i) 


loo2 |< <_< (S.87 


ao.@ [ise jaca 


19 





where 6 = arc cos —" ewe 18! Geraci) 


2¥-020 
O z 
and Gog) Sa ie (5.2399) 


Equation (5.37) is an exceptionally good approximation for 


the Ultraspherical polynomial in the domain: 


ele a 


2¥-A20 


Al | CS) 








The resulting condition (5.37) has exhibited, in a large 
number of tests, excellent agreement with experimental 
results. 

The stability conditions for the third order 
transfer function in two-dimensions are summarized as follows: 
1 
flez,5,Z2,) = = 

L = 1021" ~ 91237 - a2027° 


2 
represents a stable algorithm, when a, > - —=t4 if and only 


it 
Gin) So = | aro | = loo. | 


Z 


when a9 < - a hes Suntec teom: that 


og | < Cee = 


or alternatively 


on) | 


sine ary [02 + In%Fag| [2 - Far 


|oo1| s i 
where 6 = arc cos (au — | and Crs. i. 
2¥-a2 0 


die) 





Ieee Fourth—-Order Transfer Function 


For g(z,) as defined in eq (5.4) and 
h€z.) = Oo122> 2 AAS 


equations (5.2), (5.3) become 


a 


. eaten C554) 
ay = = => 
L- Qy9Zp° =-% 1 22° — G29Z7~ — 9227 


(oe) 


_ Z F Ieee pal: 
= [oo 122° + 92247] P [e="ey027 — OyoZT ed) P (5.42) 


Equation (5.41) is realized recursively in direct forms 3 
and 4, as shown in Fig 5.3. A hybrid structure realizing 
eq (5.42) is shown in Fig 5.4. 

a. Unit Sample Response 

The unit sample response 1s computed in Appendix 

I for the fourth-order filter in two dimensions using the 
Same approach as in paragraph 2a. The result is listed in 
mapees 5.3 and 5.4. 


Dae Neon OLeAbsolutre Entries of 
Unit Sample Response A 


The sum of all absolute entries of Tables 5.3 
and 5.4 will be obtained using Cauchy's first and second 


method, 1.e., 


Taal 


=~ 7 
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mig o.o: DIRECT FORM REALIZATION OF FOURTH- 


ORDER (UNCOUPLED) FILTER IN TWO- 
DIMENSIONS. 
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> = h( k,n, ) 


[k/2] K-24 


y (k=i1)! a> Ore 
on! CSE 
i=0 
ic 25) : : 
(k-itl1)! ah ee Oop Qo] 
The kei)! 11 
1=0 
e727] oe eee? 
(k-i+2) ! ce O12 Gor 4 (k-i+1)! os a O20 9 2 
tad >4 et Ite Ck—2i 1 at 
i=0 
[k/2] Se G2 
~ ; . —2m 
yo ce oo (j-mtk-i)! ai sae 
i! (k-2i)! ae m! (j-2m)! 
1=0 m=0 





TABLE 5.4: ny, = k ENTRIES OF UNIT SAMPLE RESPONSE 
CORRESPONDING TO FOURTH ORDER UNCOUPLED 
TRANSFER FUNCTION 
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Summation by columns: 








oo oo [n,/2] ni-21 2061S [n2/2] No-2m m 
aS 2 Q) Q20 (n)tno-m-1) Td 1 O02 
r i! (nj-21)! m! (n,-2m)! 
| (5.44) 
n2=0 n,=0 1=0 n= 


Summation by rows, respectively: 





oo co  j[n2/2] no-2m m_ {[n,/2] ny-21 21 
eS = Co Ag 2 (nytno-me1) ! 16 Q20 
7 m!(n2-2m)! yy 1! (nj-21)! 
; Coto 
n,=0 no=0 m=0 1=0 


(1) Necessary and Sufficient Stability Conditions 
It is proven in Appendix I that the fourth order 


Mmeoupted transfer function is stable for Gay > 0 and doz > OQ , 


m.e., 
5 ah 
o0 oo a iriy 2 fecy [n2/2] m 

> = Wi oO ee O10 (nytne-mi)! / do2 

Z Joo] [cto Gay aa m!(n»-2m) ! we 

i Qo1 

; ; 5.46) 
n2=0 n,=0 ILC m=Q 


Memrinite, if and only if 


Qeo t+ Ao2 < l - laro | - lao a | (5.47) 


ute) 





2 
Qo} 
uy 
of columns S° is written in the following form 





Mmiethe Case where 0 > go> > - and d29 > 0 the summation 


[n,/2] ny-221 a 


Cc = > ) |. Q10 Q2 0 (n,-1)! ae 
met (n,-21)! 


fio = 0 ii 
[n2/2] _ 
AL (-1)™ Snitne-mei)! O01 ae 
(nj-2) ! m! (no-2m) ! EDS 
m=0 


But the term in the square bracket defines the Ultra- 


spherical polynomial of the Cent) order, 1.e., 





ee 2/2) (pm (nitne-mi)! O01 oa (ni~itl) Oo 1 
Gi) ! » m! (no-2m)! 1 1G) ) 
(5.48) 


m=0 


and 


i y ny-21 
= Chic ee Ze - (ni-i)! ( ae ain Oo 1 
aye ee. ne TY 


no=0 n,=0 1=0 


It is proven in Appendix D that the Ultraspherical 


polynomial as defined by eq (5.48) is positive for 


e4 > 1 (5.49) 
2V¥-o2 


127 





and positive/negative for 


— = al (5250) 
2¥= 002 


for even/odd n2 respectively. 


Following the same reasoning as in the previous section, it 


ean be shown that for G25 > 0 , 





(ni-itl) e 
EC ( a > 0 (5.51) 
a2 2¥-Q92 
whenever 
ae 
O > Qo2 a Th Coms2e 


Similarly it can be shown that the summation of all entries 


of <A> by rows, i.e., 


co co (n2/2] No-2m om [n/2] mi-21i 0 
of a Qo Qo 2 eo (nitno-m-1)! Qio 29 
m! (n2-2m)! It (ny-23)! 
Mmr=0 N.=0 m=0 TED 


oo co [n2/2] ne-2m mm ni 
et a Qo Qo2 (no-m)! C(¥-a29 ) a al om 
m! (n2-2m)! le 2V=a20 
ni=0 n2=0 m=0 Go5 52) 
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and that 


(no-mt+1) 
C . > 0 
2¥-Q2 9 


ror 


Q 
— 
Oo Nn 





OR 020. = 


maesummary, for the conditions 


Q 
mn 





eae anGs — UR wCig2 2 o> a (5.54) 
a2 
Shoe ee) Grivel 0 Sea fee) ele a Ca Sey), 


which are shown in Fig 5.5, the column series Spasms tmed 


by eq (5.44) can be written in the form 


fore) oO [n1/2] : 
deo \~ 
mo) >) 2D | sion 
10 
al: Joos] 1! (ni-21)! 
mM2=-0 ni=0 
[n2/2] oe 
(nitn2-m1)! Qo 2 
SS enna S58) 
01 
=0 


which is identical to eq (5.46) and for which the conditions 


for convergence are as stated in eq (5.47). 


Leds 


Geo urivedent ostaprlity Conditions 


For COnd1tiensc: 











2 2 
0>a20>- 2% and 0 > aoe > - =H (5.57) 
and for 
2 Die 
0 = = =i and/or Qo2 < = =< (ey eye), 


which are also shown in Fig 5.5, the triangle inequality 


is applied to 


o ini /2) mi-2i 1 [n2/2] n2-2m m 
10 O20 (nitne-m1)! Ao, Qo 2 
Tne me Tn! sa ' =~ COS 
fMe=0 ny=0 1=0 m=O 
such that 


lees ly 


oe) 
wa Jerr o[™2 (  tn2/2) cn sane-mei) ! foo] (= 


016 O10 
ay » > i! ( -2i)! d m! (n,-2m)! 


ee 1-0 170 G5n5o 





It 1S easily observed from the derivation leading to 
necessary and sufficient stability conditions, that the 


conditions for convergence of eq (5.593) must be 


[aro | + | ooo | <l1- laio| o |aoil. (5.60) 
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Eq (5.54) 


ey (5.55) 





(5,57) 





Eq (5.58) 





Fig. 5.5: GRAPHICAL INTERPRETATION OF 
INBQUAEIMIES (5:54, 55, 57, 
58). 


phew 





fee COUPLED FILTER 
imeeecccond Order Coupled (Bilinear) Filter 
The transfer function of the second-order coupled 
filter, which is also referred to as bilinear case in 


Refs [13] and [24] has the forn: 


mi 21522) _ a C561) 
l= 01027? - 009122? - 01127723" 


Equation (5.61) can be realized in direct form as shown in 
Fig 5.6. It is noted that for a unity numerator the 
realization in Fig 5.6 represents all direct, as well as 
the canonic forms. 
a. Unit Sample Response A 

The unit sample response A , which is defined 
in Section D of the previous chapter, is computed using 
Theorem 4.1. The resulting matrix is listed in Table 5.5, 
in which the pee entry of column 1, i.e., 


Gia) oes Qo), + k ie O11 


meerewritten as 





Ji 
i aa aa d he k+1-m, 
(l-m)! m! d G10 a 


m=0 
ae ae hee 
where eC ae ( ) Seanac Lor Kk derivative ore 3 
d 1130 


We 





qela 


Using the same notation, the } entry of the a row can 


be written as 
k 
oO a ‘i Ot ak 


Ck-m)t mt aot) 


CS 62) 
d G01 


b. Sum of Absolute Entries of A 


The sum of absolute entries of A is equal to 


1G yy) pel m 
2 ghz nitn2-m 


Q10 Q13 
5 = rea - Qo) Coe ce) 
y y ) ee act 


ne=0 n,=0 





To identify convergency conditions of eq (5.63), both the 
conditions for convergence are found for each particular ng2 
o@oeO < ny} < @ , and the sum ae to which it converges, as 


well as the conditions for convergence of 
oc 


os 
Ne 


Nn =0 


This approach follows Cauchy's second method of matrix 
Summation. In particular, the sum of all entries in column 0 
of {h(n;,n2)} is shown in the corresponding section C of 


feeendix I to result for ai, > 0 in 


Jue: 





which converges if |aio0| < 1 to 


Il 


oh ae (Le eee, 


Si, the sum of entries in column 1, 1s proven in Appendix K 


memequal, for |o%,,[| <1, to 


S, > ih Q190%01 + O21 
— l- aQi9 Oso 


Following the same procedure, the S. becomes 


j 
a ( Sag.co1 On) 
35 i= Q10 re 0 aon | CORON 





But the sum of all absolute entries of A is defined by 


No =0 


which is now identified for Lore | < 1 to equal 


oo 
, 7 
+ 
Tat eae ors > (Saeger a1 + go, ) (5.65) 
no =0 
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But the right hand side of eq (5.65) converges if and only if 


oY 





S EaCTTTae 2 51 Ccnice) 
= 110 


me Discussien 
It will be shown in the next chapter that another 
necessary stability condition for the bilinear case is as 


follows: 


Con = Oe 


ere: < l Co-67) 








Moreover it is known [13] that the conditions in eqs (5.66) 


ema (5.67) are also sufficient. 


LSS 





) 


Loe: 





Fig. 5.6: TWO DIMENSIONAL SECOND ORDER 
COUPLED FILTER IN DIRECT FORM 
REALIZATION. 
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Vw ieGeoonny son SOUPPICEENT STABILITY CONDITIONS 
IN N-DIMENSIONS 


A. INTRODUCTION 

It was shown in Chapter V how conditions for absolute 
convergence of all unit sample response entries, which were 
shown in Appendix C to be identical to BIBO stability 
conditions, are derived for the third and fourth order 
uncoupled and second order coupled case. The removal of 
Mmeolute signs enclosing h(n,,n,) in the Taylor expansion 
of a two-dimensional transfer function is possible in general 
by application of the triangle inequality or whenever all 
coefficients in the denominator polynomial are positive. 


This will be shown in the following paragraphs. 


im Necessary and Sufficient Stability 


GConaltvonse fOrl-osauive Characteristic 
Se i De ee ann Siemainennnn 
Equation Coot uenetents 19) 17e- 1 mensilons 


The summation of <A> entries is easily made for 
the case where all coefficients of the two-dimensional 
characteristic equation are positive. In fact, the sum of 


all entries is obtained by observing that 


ies) 522) = C615) 


M 
-Mm, -m 
L = 4 2} AN 2 
Mm )MN2 


m,-0 mz —0 


(m,t+m> A 0) 
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M, M, Q 


-m,; _-m 
fez, 5Z5,) = iy » » ne Zz (652) 
Q= m,;=0 m,=0 


(m,+m,#0) 


The Taylor expansion of H(z,,z,) is defined in eq (4.51) as 


pe 


» Cin cine) aqme aaa (6.3) 


0 no=0 


oo 


Ny 


“A 
a 


Therefore, each and every h(n,;,;n.) can be identified from 


eq (6.2) by computing the binomial expansion and combining 


terms of equal powers in z, and z,. Moreover if all a > 0 


Mm Mo 


then by theorem 4.4 h(n, ,n,) > 0 and it follows directly 


that 


q , Jh(ny,n,)] = » » inne min Coa) 


n,=0 n,=0 n,=0 n,=0 


fests proven in Appendix C that H(z,,z,) is BIBO stable if 


and only if 


» % UstiGewuryey yy ican C655) 


n,=0 n,=0 


ILS 


or by application of eq (6.4) if and only if 


J yy h(n,,n2) < » Coe 6) 
N} 


=0 n2=0 


But the sum of all unit sample response entries is also 
identified from the right hand side of eq (6.3) for 


m= Z>5 = 1 , which by equating eqs (6.2) and (6.3) equals 


oo 


y » h(nj,n,) = » , y —e Ca) 
,=-9 mM, 


m—-U0 n -O0 g=0 


ees Ore it 1s stated that H(z,,z,) is BIBO stable for all 


On | m, 2 ieee and enly af 
— My Mo 


| Q 
q . > Om,m, < co (6.8) 


m,=0 m,=0 
(m,+m,#0) 


— 
it 
© 


But necessary and sufficient conditions for convergence of 


eq (6.8) are easily identified as 


M MM 
) yl Omm, < 1 (6.9) 


m,=0 m,=0 


(m,+m, #0) 
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ewok lorent Stability Conditions in Two-Dimensions 


The sum of all <A> entries is by the triangle 
inequality less than or equal to the sum of all entries of 
the unit sample response matrix in which each coefficient 


is replaced by its absolute value, l.e., 


co 


roe) co My Mo 
L 
y q Jh(ny,n,)| < yy » » Ca | CEO) 


n,=0 n,=0 L=0 m,=0 m,=0 
(m,+m,#0) 


But the right hand side of eq (6.10) converges if and only 


ibq 


My Mo 
q a lorem, . : ? 


m,;=0 m,=0 


which by eq (6.10) is also a sufficient stability condition 


Meret Z,,Z.). 


3. Necessary Stability Conditions in Two-Dimensions 
By triangle inequality 


co co 


y » Jh¢n,,n,) | iz | » » IgG okey ote, (6.12) 


n,=0 n,=0 n,=0 n»,=0 


which, by eq (6.7) is equal to 
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fore) co co M, Mo 


x 
|) J h(ny,n2)| = |) % * a (om) 
10955 (IQ ale} 


ni3=0 n,=0 2=0 m,;=0O m2=0 





Phe right hand side of eq (6.13) converges if and only if 


My Mo 
> ds 
myMmo 


m,=0 m2=0 


a (6.15) 





Equation (6.14) represents necessary conditions for stability 
me( 21,22). 

For the following sections and in Chapter VII, a 
mathematical generalization of the above conditions to 
N-dimensions, improved necessary conditions and a method to 
determine necessary and sufficient stability conditions in 
N-dimensions are presented using the theory of modern 
algebra as in Ref [53]. Notice from the unit sample response 
corresponding to the third order uncoupled transfer function, 
that the signs of a9, OY a;)9 have no influence on the 
convergency behavior of <A>. This observation can be 
mathematically formulated in general as "symmetry condition" 


in the N-variable coefficient space. 


See COEFFICIENT SPACE 
The expression (S1,-+-+58)) denotes, using the notation 


established in Chapter IV, a fixed k-tuple of vector 
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indices, each of dimension N, i.e., for the third order 


uncoupled two-dimensional transfer function: 
Cie cmo) e- C lop (Om; (2,0) } 


The characteristic equation Q(z), of an N-dimensional transfer 


function L(z), where 


S. 
OGZ) = 1 » A. z °: 8 2S), 


is uniquely defined by the arbitrary k-tuple of real Nunbenet 
(Ai,...,A,) and by the ordered k-tuple (Si 526+ 55)) of vector 
imaices. 

For a fixed k-tuple of indices (Si5..+55)) ie eer oc cull 
k-tuples (Arye. + A.) determines the coefficient space for the 
polynomial Q(z). Each k-tuple or 'point' in the 
k-dimensional Euclidian space oe is therefore associated 
With the polynomial of the form of eq (6.15). 

it 1S said that the point (Ais++-sA,) is a ‘stable point' 
if Q(z) associated with the k-tuple of real numbers corre- 
sponds to a stable transfer function. 

The set of all stable points in ites is called the 
Bremen of stability!’ in IR* determined by the ordered 


k-tuple (Si5+-+5S)) of vector indices. 


1430 


Example ole 


The stability conditions for the second order coupled 
two-dimensional transfer function derived in the previous 
chapter define the region of stability in the three- 
dimensional coefficient space (A,,A,,A3;) as shown in 


fmmes G.1, 6.2, and 6.3. 


meee NECESSARY CONDITIONS 

In this paragraph, two necessary stability conditions 
in N-dimensions are derived, namely, that the absolute sum 
Beeche leading coefficients and of all coefficients in Q(z) 
must be strictly less than one. The term "leading 


coefficient" is defined in Appendix A. 


Theorem 6.1: Let (Ayye++ AL) be a point in the coeffircient 
Space and let A be the set of A, which are leading coeffi- 
cients of the polynomial associated with CArs++e Ap). Ef 


(Ay,+++,A,) is stable, then 
pa 
it 

A-eA 


Proof: See Appendix M. 


<i 





Example 6 wk: 


The two-dimensional full quadratic transfer function 


has the form 


ee SS * (6.16) 


2 
L = O30Z3 — Go1Z2 — 4112122 — G20Z1 — Ao2Z2 


144 





The set of leading coefficients is identified as 


Pearce) (2.0). (0e2) | 


A necessary condition for stability of eq (6.16) requires 
that the absolute sum of the leading coefficients, which 


are 411, G29, and G2, 1s less than one, i1.e., 
lora + Geo + Go2| < 1 
Another important result regarding necessary conditions for 


stability of an N-dimensional transfer function is: 


Theorem 6.2: Let (Pise.+sPy) be any N-tuple such that 


a = 0,1, or -1 fori:=1,...,N. Define the k-tuple 


(Pis+++sPy) by 


D. = (riy... yr) h Cena) 


rh 
itien, if (Arse sAL) is a stable point, it must be true that 


BEE pon 4 ey <7 dl (6.08) 


Proof: See Appendix M. 


145 





Example Go: 


For the full quadratic transfer function of the previous 
example, let (r,;,r,) = (1,-1). Then the 5-tuple 


Pa >P2 >P3>P, »Ps) is defined by 


oy 2 Chey SEO) San Sa 
pal (ge) Ue ce eae 
Dg 2 (geal) mene seal 
Py = Cujcny ow 1 
Ps = ean Se? Hf 


Therefore, by Theorem 6.2, if 
(Ar,.+-,A,) = i a 
is a stable point, then 
Gages fo1e> 8117 * G20 * Goo < 1 


Whenever all rn. = iImMion 1e—si....N the following corollary 


can be stated: 


Corollary i: Tet: CAr,..+,A,) 1s a stable point, then 


JV een Lar <~ i OSVRAL S|), 
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Theorem 6.3: Let (Pis.+- Py) be any N-tuple such that 


a > mlseer -l fori = 1,...,N. Define the k-tuple 


(Pis+++>sP,) by 


Se 
al 


Pi = CS er ave, C620) 


Miven, if (Ars.+.9A,) is a stable point, it must be true that 
|piA 7 ae: PA, | soil C6221) 


Beeot: See Appendix M. 


Finally, following the above reasoning, corollary 6.2 is 


formulated as follows: 


Sevollary 6.2: if CA1s.++sA,) is a stable point, then 


Us BF one ae < tL (6.22) 


x | 


Bee oUFFICIENT CONDITIONS 


Theorem 6.4: Let (Ars++-,A,) be a point such that 
k 
> tas eal Ce 
aja 


then (Ar ,.+.,A)) is stable. 


Proof: Follows directly by generalizing the derivation of 


section A2 to N-dimensions. 
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Example Gee Ue: 


Pesurfacient condition for stability of L(z,,z,) as 


defined in eq (6.16) is: 


fonotet lool + ergot leo) + |ago| < 1 


A summarizing discussion of above stability conditions is 
presented in the next chapter. Note that necessary stability 
conditions and sufficient stability conditions are listed 

in tables 6.1 and 6.2 for two-dimensional filters up to 

Sixth order and for selected N-dimensional filters, 


respectively. 


meee YMMETRY CONDITIONS 


Theorem 6.5: Let 6 = (6,,...,6,,) be any complex vector such 
that 
1) 15, | = 1 for alli Co.24) 


meeeecor i= 1,...;K the numbers Ps defined by 
Si 
ee == (S15-+-55y) (6.25) 


are all real. (6.25) 


Then a point (Ary++-sA,) 1s stable if and only if the point 


(Als. ++ Ay) is stable, where 
aE jo\e for 4) S doggone (6.26) 
Proof: See Appendix M. 
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Example 6.4: 


The third order uncoupled two-dimensional transfer 


function 
a ae 
1 = 41021 - G01Z22 —- 42021 
mao = —L 
Pi = (-1,-1) S129? = =] 
Bere Caan) 82? = ag 
55 = Cone Soo 


Thus, if (A1,A2,A3) 1s Stable, then the necessary and 
sufficient condition of Theorem 6.4 prescribes that 
(-A,,-A2,A43) is stable. This implies that the stability 
region of L(z1,2Z2) 1S symmetric with respect to the (A)) 
and (A2) axis, which was already observed in Section A of 


this chapter. 
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1 
Transfer H(Z1,Z2) = = a me 
Function 1 - ©1902) = 09322 =~ Go9Zj 









STABILITY CONDITIONS 







Necessary looo] <1 Theorem 6.1 








loro] + [aor] + deo <1 Theorem 6.2 
loro + Gor + Goo] <1 
|-c19 + Gor + Gao| < 1 
Theorem 6.3 
loro — Oo1 + Goo| <1 
|—o10 — Oo1 + Gao| < 1 
sufficient loro} + |aoail+ loso] <1 Theorem 6.4 


Table 6.1A: Sufficient and necessary stability conditions 
for the third order (uncoupled) two-dimensional 
beats rene | uncon . 


NG 





Transfer H(zi,Z2) = 
Function 


00 P= 
nA 2 


1 - Omnm> Za 
m1=0 m2=0 
(O<mi+m, < 3) 





STABILITY CONDITIONS 


Necessary loso + G21 + O12 + Oo3| <1 Theorem 
6.1 
[oro + Gor + G11 + Goo + Go2 + Gao + O12 + Azo + Oo 3| 
ae 


|—cir9 + 92 + G11 — Geo + Ao2 + G21 + Giz — Azo + Oo 3 | 


=k Theorem 
6.3 


lo10 — Go1 — G11 + Gao + G2 — G21 + Ai2 + O30 - Oo 3| 
ol 
|—c19 — Gor + O11 + Geo + Ao2 — O21 — G12 — A30 — Oo 3| 


a 


3 3 
Sufficient , 


m ,=0 M.=0 
(O<mi+m2 < 3) 


Theorem 


6.4 


loamy 





Table 6.1B: Sufficient and necessary stability conditions 
for sixth order (coupled) two-dimensional 
iecwiomem= mune tl On. 
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Mieanstfer Furrction 
F (2,22) = 2. 








re . on 
“ 4 - ©, = Sere 
m, 
od « 







SABIE eC ONPETTONS 






Necessary: 


Theorem 6.1 






| “iy, ah + oN yas = deat Cie dt di) mE do, 1 < 4 


i. ' ; : 5 } 7 s) 
Ang + Aer, a cide ao en sie 1 Vins 7 Ae a) dey [4 a 


QoL 


‘ 09 Sr eo, + Car Ore af di, i Py [<3 


os ings + lees of #hscy 4- ai, = cr - uo Ec ie 
Theorem 6.3 


cs Ae aE oe + Avert + dy, mr aie - de, Ne as 


] ' 
+ Nios + ace + Boy - & Ip & Slo ~ Apy }Z ne 


+ Sao + Sort ae = dn, a Chat oF ode \ 





sufficient: 
2, i ae 
ee 


Theorem 6.4 
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Table 6.2: Sufficient and necessary stability conditions 
for full quadratic three dimensional transfer 
funewron. | 
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Viti li cusmakt oANoeo Ure tC LENT CONDITIONS 
IN N-DIMENSIONS 


A. INTRODUCTION 

Pestability region 1s said to be convex if the loci of 
all points Reece two stable points lies inside the 
stable region. It is well known that the stability area 
for a second order one-dimensional transfer function satisfies 
convexity, but it can be ‘shown that the stability region of 
the third order one-dimensional transfer function is not 
convex. It is, therefore, conjectured that the bounding 
surfaces of a multi-dimensional higher order filter stability 
region cannot be described in terms of plane surfaces which 
explains why the following method to derive necessary and 
Seerrcient Stability conditions is rather complicated to 
apply. The advantage of this approach lies in the fact that 
the method must be applied only once for a particular type 
of transfer function, which results in a set of stability 
conditions. These conditions can be easily applied by the 
design engineer. 
meets tHOD TO DERIVE NECESSARY AND 

MUEFICLENT STABILITY CONDITIONS 

The method to derive necessary and sufficient stability 
conditions consists of two steps: 

1. Derive necessary and sufficient conditions for one 


variable polynomial with complex coefficients using 
Rouche's theorem. 
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2. Apply the following theorem and corollary to state 
necessary and sufficient conditions for N-dimensions. 


Theorem 7: Let (S150 5p) be any N-tuple of complex 


mimbers such that r = max 1S | <1. Then Q(21,..-,2y) 
@mstaeble if and only if for all 6 the polynomial 
Bn oN 
qs 62) = OiGZ a Z =e C7 e2) 


1s stable. 
Proof: See Appendix N. 


Corollary 7: The polynomial Q(z) corresponds to a stable 
transfer function if and only if for all complex numbers 6 


such that 
max |S. | eal qs 2) = Q(261,..- 526) Cia) 


memsotable. 


Proof: See Appendix N. 


Smee xAMPLE 7.1: FIRST ORDER TRANSFER FUNCTION 
IN N-DIMENSTIONS 


Using the notation of the previous chapter the first 
order case in N-dimensions is defined as follows: 


N 


Ss. 
Oe) a= see y A. z 2 Gig) 


1=1 


| 


where 


ie Cae. ..0) 
ees (0 eo) 
ae: *GeMo.....21) 


N 
q(z) = 1 - » Peon 2 commmmtax |S) 0-2 


<i Gra) 





Mmiemby corollary 7 Q(z) is stable if and only if eq (7.4) 
holds for all 6 , where max |6;| = 1. For 6, = #1 


eq (7.4) results in 


[Ai| + ... + [A oe! Gi) 


| 


Which is a necessary condition. But by Theorem 6.4, eq (7.5) 


158 





is also a sufficient condition. In summary the first order 
N-dimensional filter is stable if and only if the sum of 
the absolute coefficients is less than one. 
D. EXAMPLE 7.2: N-DIMENSIONAL TRANSFER FUNCTION WITH 
Bost itive COBFFACIENE DENOMLNATOR POLYNOMIAL 


men O(Z) is defined as: 


and all as > 0 , then a necessary and sufficient condition 
ah 


for stability is that the sum of all coefficients is less 
iran one . 

This follows immediately by applying Theorem 6.2 for 
a = 1 to derive the necessary condition and Theorem 6.4 
bor the sufficient conditions, which for all coefficients 
positive are identical. Note that the stability condition 


for the above cases as well as for the bilinear case 


(without proof) are tabulated in the following Table.(Table 7). 


EB. SUMMARY 

several powerful methods have been developed to solve 
the stability question of N-dimensional recursive digital 
filter. Ultraspherical polynomial and derivative operator 
methods were used to derive stability conditions for specific 


low-order two-dimensional cases. Methods for the derivation 


JUD) 


of sufficient and necessary stability conditions in two- 
dimensions and in N-dimensions and for necessary and suffi- 
cient stability conditions in N-dimensions were presented 
and applied to compute low order two-, and N-dimensional 
stability conditions. 

A more complete derivation of stability conditions for 
low order two- and N-dimensional recursive digital filter 
and an arrangement in tabular form are left for future 
research. With the metheds presented to solve the stability 
question the second important part, l1.e., approximation 
of a given specification by a transfer function, will be 


investigated in the following chapters. 
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foeel., PROPAGATION OF TRANSFER FUNCTEON COEFFICIENTS IN 
fie UN aoe LE RESPONSE OF N=DIMENSTONALDSEILIER 


A. INTRODUCTION 

In the previous chapters, stability criteria for 
N-dimensional filters and implementation forms were developed 
which will be applied in the following chapters to derive 
time domain design techniques of N-dimensional recursive 
Seeital filter. The design of a recursive digital 
filter is achieved in two steps, (1) choose filter coeffi- 
cients to approximate a given response, and (2) ensure that 
m@eeresulting filter is stable. In Chapters IV, V, VI, and 
VII, several techniques to solve the stability question were 
developed and stability criteria for N-dimensional and 
specific two-dimensional recursive filter derived. In the 
following chapter the nonrecursive combinatorial formula of 
Theorem 4.4 will be utilized to relate transfer function 
coefficients with the entries of the unit sample response 
in N-dimensions. iat will be specifically shown that :the 
coefficient of the characteristic equation acy) appears 
only once in linear form and with coefficient one in the 
coefficient h(n) of the corresponding unit sample response. 
It also appears in nonlinear form, with non-unity coeffi- 
Gents in all h(j), where n<j. 

This important correspondence will be used in Chapter IX 


to construct an algorithm for the extraction of an all-pole 
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transfer function from a given unit sample response. The 
extraction of a rational N-dimensional transfer function will 
be achieved in the same chapter utilizing a method which is 
essentially a Padé approximant extended to N-variables. 
Both of these methods will be incorporated in an algorithm 
extracting a parallel arrangement of low order sections. 

In Chapter X the time domain design step 2 will be 
achieved by approximating a given unit sample response 
Magen) } , where O<n; <7, O< i<N by fh(n)} , 


Omere QO<n. < k and then using transfer function 


ees? 
extraction methods. In step 2 of the design procedure, the 
approximating transfer function is tested for stability. The 
Quality of the time domain design is evaluated comparing 


the spectrum of {h(n)} to the spectrum of the extracted 


mrmansiter function. 


B. PROPAGATION RULES IN N-DIMENSIONS 
It has been stated in Theorem 4.5 that the Taylor 


expansion of 


H@) = ———_—_, 


Le Ota — 5 


where 


Oe = 5 » a(T) z + 


w= 


ae 
(0<i}<@) 
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has the form: 


The coefficients h(n) are related explicitly to the aT) 


by theorems 4.4 and 4.5 as follows: 


— + + = =~ 
¥ \ cc CommNGnicic oa 
h(n) = CS) 2 ) 
ey, : : Sy : 
t= t= 
So : Sy 
(8.1) 
where 
te, 80 2 ea Sere age = (=) Come) 


In the next paragraphs, three cases for the appearance of 

a particular coefficient as y » where (0 <p << Kk) ,; 

in the coefficient umecetae h(n) of the Taylor expansion 
Semiicz) are investigated, namely (S..) Pai ars (S.) = (ns 
and (S,,) Git) 
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1. Propagation Rule for az) » where CS a) 
D p 





Tac noe tation (S_,) SEChvmes detuned as 


(Ss. > (n) 


or (S.? > (n |) 


Py N 


For example, for (S (2,1) the inequality 


LS) 
Pp} D2 
(2,1) > (ni,;n2) designates the shaded area in the 


following graph: 





ma eae OD C0, (ge) Co soe) 
Zia a va A ve 
1 Claret ma@ine) Cie 3)- (1.4) Ux.. 
ve 7 a Ze 
2 Oe C2. 2 2) 02.3) (2.4) 
ee 
a Con ebrGent), (352)5 (33a 3,4) 





G5 


Propagation Rule I: The coefficient acs) of the N-variable 
_ see os 

characteristic equation of H(z) does not appear in h(n), 

whenever Ce - Ch) ae ie remenGn)) 1S iene aoe coefficient 


Siethe Taylor expansion of H(z). 


Example 8.1: 
ror 
dl 
ECZ5 5 25.) =~ 
—1 = “= —~s_-4 
ISO An = o122 i 3021 a 5421 22 
C836) 
the corresponding unit sample response {h(n;,n2)} is 
computed from eq (4.21). In the case where (n,,n2) < (5,4) 


eq (8.2), which becomes for the above transfer function 


toy,0) 6420 * tg 4) 75D + tog 9°30) + tog yy C524) = (ni yng) 


can only be solved, if and only if, te. oy = 0 
3 


2. Propagation Rules for ag) > where (S.) = Ga) 
I a) ee 


The Netatven (S) = (nm) is defined as 


oo? = (n)) 
and cS = iC), 
P2 
and ad = (ny? 
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if (S,,) ei. | them there exists one term in h(n) where 


te equals one and all other t's are identically zero, 


therefore eq (8.1) can be rewritten as follows: 


= Sp Dp terms of 
mem) = el + a 

oa) 
te 0, 0 ( D order < Sp 

P te =a) 
p 

= terms of 

oS)” = 

P order < 5S 


In this case acs ) appears in linear form with coefficient 


ene. lhis observation is stated as follows: 


Eeeparation Rule JI: The coefficient a5 ) of the N-variable 
Pp 


characteristic equation of H(z) , appears in h(n) in linear 


form with coefficient one, whenever (S$) = in) 
Example 8.2: 
The coefficient Ars yy = O54 of the same transfer 
> 


function aS in example 8.1 appears in h(5,4) in linear 


form with coefficient one, since 


(1,0) +t (0,1) +t (3,0) +t Gate ase) 


a0) (0,1) (3,0) (5,4) 


whenever tc L) = 1 and all other t's are zero. The terms 
> 


of order < (5,4) are identified as: 


ey 





umes «60,1 ~~? FS 3,0) > * Fes =O 3 
and 
mies” °C 00,1) > * =? «683,09 79 § eS uy = 9 
Therefore, using eq (4.48a) 

he ,4) = as, + 105 an oi G39 + 112 O16 Oo 


If the i diagonal cut of the unit sample response is 
feeamed to be the set of indices (n) = (Niy++.sMy) , such 
that 

N 

> n, = I . (8.7) 

1=1 


then the following corollary to Propagation Rule II can be 


stated: 


maepesation Rule IIA: Consider the first diagonal cut of 
{h(n)} , then each entry corresponds directly to a first 


Sreaer coefficient are) of the N-dimensional characteristic 


Samation of H(z) . In fact, for 
N 
e=al 


h(n) = as ) » whenever (S.) ar) 
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Example Saye 


The unit sample response entries h(n) along the first 


diagonal cut corresponding to the transfer function 


H(z, BZ) = 


a 


-1 -1 -} -) -] -5 -4 -3 
1-100Z21 — 91022 - Q001Z23 — 11021 22° — 54321 Z2 23 


are Mga, 0) = i920 
heGOe , 0) = @9010 


IgiCiel (Ope => 9001 


This follows directly from eq (8.2), which can be written 


as follows 


fame) 22029? * tq 49802220) * Teg 91602052? 


a ¢1.1,0)'t242 + t¢5 4,3) Oe? 


=e GlesO- 0) GO. 0) 1s C001) } 


(End of example) 


appears 


with coefficient one, whenever (n) = WS, , then 
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g Sp S ELITE Sse 
h(n) = a + == 
S < S 
te 20 ( >? order < W ( a 
p = 
te W 
p 
 . W re terms of 
: (S, order < W (S) 


Therefore, the following rule can be stated: 


Mm@eeagation Rule IIIT: The coefficient as) of the 
Pp 


N-variable characteristic equation of H(z) appears in h(n) 
h 


to the we power with coefficient one, whenever (n) = We 
Example 8.4: 
The coefficient a(5,4) = as, of previously used 


transfer function appear in h(10,8) with coefficient one 


to the second power, which can be identified from 


PG 0) nF ey Oem + t( 3,9) 6320) + Si) Sb =. C0, 89 


where for tes Lb) = 2 all t's must be zero. The terms of 
> 


order < (10,8) are computed from: 
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meme | Coe) ~ © 3 *¢3,0) ~ & ~~ ¥esjuy = © 
mang) ~ / aie) = Seno) 1 om C5 gu) ee 
Meme C1) ~ & 30 8,0) = * = e5,ny = 2 
Meme > * 8 (0,1) ~ * 43,0) >> £tes,u) = 2 
Memo) Ci‘) ~hCfCUCUC8,0) ~ ® 385 uy = 2 
Thus 
2 18! 10 8&8 16! an 
mol0,8) = as, + Torey 10001 F ayer %10%01030 


14! h 8 2 12! 8 3 
+ Brut ar 1000 1030 i Br ar %10%01020 


gt 2 icp ae: 
+ FT yr A100 1030054 + Spay 410001054 


h 


meopagation Rule IV: In the me diagonal cut of mem)! 3 


each entry contains an yeh 


order coefficient aca) of the 
N-dimensional characteristic equation of H(z). In fact, 


min) = = + (terms of order < Sy 7 where (3, = (n). 


lye 


3. Propagation Rules for a5) » where (S)) << Sy 


p 


The notation (S.) “Gn eo sedetined (as an ireorems td) 


as 
(S, > Gn): 
cameo...) < (lo) 
La 
andi (S. -)u< (ny)? ; where the case 
Py a 


(S,) =e eomerpileidl: sexceluded: 


For example, for re oe C2) tie aeamMe dial ty 


(2,1) < (nmi,;n2) designates the area, which is shaded in the 


following graph: 





eae OMC OO meu. dt) C052) (0,3) -CO,4) 


er) tender Cl Zeel. syle) 





ie 





Example co One 


The set of coefficients {acs 3 » where i1,tli, = 2 
1 


dey) 


appear in the coefficients of the Taylor expansion of 


L 
eZ 5Z2) = 
1 -1l1 —]1 2 -—2 -5 —4 


all = = 
Peco 21). Ogi2Z2 — @)12) Zo — 202), = 9222 = A542) Zo 


(8.8) 


Mmene the second diagonal cut, i.e., in {h(€n,,n,.)} , such 


that n, +n. =2 =, #42Since 

£67 96229) * tg 3) 6022) Ftp 4) C1od) t teo gy 6220) t tg 9) (052) 
* tes yy So) = (2,0) 

t (1 9) (120) * tog 4) (0D toy yb # toy 9) (2:0) # tg 96052) 
* £5 yy 654) SiC) 

t (19) (120) # tog 4) (0D) tte yy Ged # toy 9) (20) * tq 9) (022) 
*t 5 yy ho) = (0,2) 


Consequently, we obtain by inspection 


vd 
m2 50) = doo + Qj,y 
imc] , 1) = 4) Tt 2020001 
m0), 2) = Qo2 + Oot 


ae 3 


3 


Eos (S.) “wo eeence srolllowing midies Can be stated. 


Propagation Eule V2" Uhe coetficient acs) of the N-variable 
= P = 
characteristic equation of H(z) appears in h(n) for 


(S,? mcr e, whenever equation (8.2) is satisfied. 


Meemasgation Rule VI: The coefficient acs} of the 
_P 
N-variable characteristic equation of H(z) appears in h(n) 


for ws) Gm winenever eqtation €8.2) 16 satistied, 


fee © ROPAGATION DIAGRAMS FOR TWO-DIMENSIONAL TRANSFER 
mee. TONS UP TO THE SIXTH ORDER COUPLED CASE 


The Taylor expansion coefficients of two-dimensional 


mee order coupled case, l.e., 


al 
H(z — ) atl = -1 -1 =2 -2 3 
I= On¢29 = Opieo —- 01121 22 =— 02021 ~- Ao2Z2 — W392) 


C822 

a ee a2 5 
— 0212) 22 ~ 41221 22 =~ Ap 3Z2 
Will be investigated in this section. It was shown in Section 
D of Chapter IV, that in two-dimensions the unit sample 
response matrix A is obtained by writing the Taylor 
expansion of H(z,,22) as matrix product in the following 
manner 


(1 zo: Zo" Zo” veeyt 


[>| 


Meza ,Z2>) = (1 2. 21° z1° ae 
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where 


IL INGOs) icin) incOines) 

mG. 0) ] og ol ke Ges 743, Gl. 3) 

A =1n(2,0) (2,1) (2,2) h(2,3) 
ines, 0) Gs) hie 3.2) nts, 3) 


If for a moment the case is considered, where in eq (8.9) 
all coefficients are zero, except 319 and Oo}, 1-€e., 
it 


I Zy 922?) SS — (8520) 


Ll = 1021 - 0122 


then the corresponding unit sample response matrix A is 
Shown in Fig(8.1). The coefficients a315 and 9; ;, which 
are by definition (see Appendix A) also the leading coeffi- 
cients of eq(8.10), appear by propagation rule IIA in h(1,0) 
and h(0,1) in linear form with coefficient one, and 
"Dropagate" by propagation rule V for all (n,,n2) > (1,0) 
and (n,,;n2) > (0,1) , respectively. These considerations 
are easily identified in Fig 8.1. 

An abstract method to graph A is shown in Fig 8. , 
in which each entry of the unit sample response matrix is 


replaced by a cross, and single lines outline the sector of 


[>I 


in which the coefficients a), anda), "propagate", 


1.€., appear whenever eq (8.2) is satisfied. 
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Using similar symbology the leading coefficient 


propagation of the third-order uncoupled transfer function, 


1.€.5 
: al 
m.Z1,22) =  ————=—._ | tana C8 adele) 
Oto 2 FO 2 eo oy 
meesnown in Fig 8.2a. The propagation of one wench 


follows from propagation rule VI is outlined and identified 
by the coefficient aQ29 raised to ae W power. 

Figure 8.2c shows the leading coefficient propagation 
for the second-order coupled (or bilinear) transfer function. 
It is observed that a propagates in accordance with 


h entry along the main 


propagation rule VI from the we 
Seasonal for all (ni,n2) > (1,1). 

Figures 8.2d and 8.3e display the leading coefficient 
propagation of the fourth-order coupled (or full quadratic) 
and the sixth order coupled transfer function, which follows 
from propagation rule IV, where I = 2, 3, respectively, 
and propagation rule V. 

The following table relates two-dimensional transfer 


function and leading coefficient propagation charts up to 


the sixth order coupled case. 


D. SUMMARY 
The six propagation rules developed in Section B from 


mee non—-recursive combinatorial formula of Theorem 4.4 
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establish a relationship between the coefficients of a 
transfer function in N-dimensions to the coefficients of 
an N-dimensional unit sample response. The rules state 
explicitly that a transfer function coefficient ary 
appears only in h(i) in linear form and with coefficient 
Miewand propagates in all h(n), where (n) > (i) . 
Conversely, if from a given unit sample response it is 
desired to identify the transfer function coefficient Q(T) > 
then by propagation rule II entry h(i) contains aT) > 
as well as nonlinear terms composed of a(3) sesuech. Eitare 
mee . But, if the unit sample response of the transfer 
function H. (21 ,Z2) having tarzy} coefficients, where 
oot , is known then aT) is equal to the difference 


nel 


ef h and the 1 entry of the unit sample response 


(i) 
corresponding to H.(Z1,Z2). Moreover, it will be shown 
in part B of the next chapter that this idea can be imple- 
mented to extract all (5) coefficients, where 

N 


1=1 


et the same time. 
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FIG 8.2b: PROPAGATION OF 
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Pema linc LUON DOr TRANSEER FUNCTIONS FROM A 
GIVEN UNIT SAMPLE RESPONSE 


me =LNTRODUCTION 

It was outlined in Section D of the previous chapter 
how the propagation rules can be used to identify a transfer 
function coefficient from a given unit sample response. 

This technique is elucidated and is used to develop an 
algorithm which extracts an all-pole transfer function from 
a unit sample response in N dimensions. 

In Section B, an extraction method for rational transfer 
functions is presented which is essentially a Padé approxi- 
mate generalized to N-dimensions. This procedure is shown 
to require the simultaneous solution of a number of linear 
difference equations, which will preferably be done using a 
general purpose digital computer. 

Finally, in Section E, an approximate partial fraction 
expansion in N-variables in terms of all-pole and/or rational 
low order transfer functions is presented. 

A measure of quality for the extraction methods will be 
established in the next chapter, in which the above tech- 
niques are incorporated to solve step 1 of the time domain 


design procedure. 
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B. EXTRACTION OF AN N-DIMENSIONAL ALL-POLE TRANSFER 
PUNE TEON FROM A GIVEN UNIT SAMPLE RESPONSE 


The given N-dimensional unit sample response {h(n)}, 
where 0 < gy aa and 0 < i< N is chosen to be normalized, 
fem, nht0) = 1. The transfer function coefficients having 
merst order, 1.e., the set of 2(T) 7 BOM Whilchmashie ms Wine oO. 


mieeeindices of (1) is identically equal to one, for example, 


{arsy} 


mee Cl0....,0)°"(0,1,0,...,0)°°°°**’0,...,0,1,0)°°(0,...,0,1)- 
are known by propagation rule II A to be identical to the 
Corresponding set of unit sample response entries h(n) , 
taken along the first diagonal cut. 

Mite entries {h(m)} , where (n) is taken along the 
second diagonal cut contains, by propagation rule IV, all 
second order coefficients, in linear form with unit coeffi- 
cients, as well as first order nonlinear terms. To identify 
mmeeset of second order coefficients, i.e., faczy te ; 
miewunit sample response entries {h;(n)} corresponding to 
fa@eetransfer function Hi(z) , which has only first order 
coefficients, is computed and the subtraction 
Min, (n)} = {h(n)} - {h,(n)} is performed for every term 
along the second diagonal cut. The resulting set of 
Seetticients Ah2(n) is identical to the corresponding 


{acsy}2 . It 1s apparent that this technique leads to an 


It 





itera 


algorithmic approach in which the k order coefficients, 


7. tacsyty are identified from the set of differences 


Ah, (nm) = h(n) - h,_,(@) 


taken along the ee 


diagonal cut. If this procedure is 
continued until {An (n) } = {(0)} , then the all-pole 
transformation H; (2) is the exact transfer function 
corresponding to the given unit sample respone {h(n)} 

It will be shown in the next section that the above 
technique can be quite readily implemented on a general 
purpose or table top computer. It is noted that the order 
of the extracted transfer function in each variable is less 
than or equal to the size of the unit sample response along 
the corresponding mesh axis. 

SeeeecOMPULER ALGORITHM FOR THE EXTRACTION OF AN 

N-DIMENSIONAL ALL-POLE TRANSFER FUNCTION FROM 

A GIVEN UNIT SAMPLE RESPONSE 

The following algorithm incorporates, in addition to the 
concept developed in Section B, a mean square error 
criterion which will be used to terminate computations. 
Example 9.1: Extraction of a two-dimensional digital 

transfer function from a given unit 
pulse response 
The N-dimensional algorithm, which is outlined in block 


diagram form in Fig 9.1, has been implemented on an HP 9830A. 


The next pages contain an example which is designed to 
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1/GIVEN: UNIT SAMPLE RESPONSE {h(n)} , 0 <n, < I, 


TERMINATION CRITERION: MEAN SQUARE ERROR & 


2{NORMALIZE {h(m} SUCH THAT h(d) = 


emrALIZE D=1,I1Ie=1, {aps} {0}, {Ahi (n)} = fh(n)} 
ror 0 ae S “= 


31 EXTRACT FROM {Ahy (n) } THE CORE UCLENIS aay = Ah, (n) 
dhe! 


ALONG THE I iain CUT ( FORALL 0 <n, < 1) 
41 CONSTRUCT N-DIMENSIONAL CHARACTERISTIC EQUATION D , i.e. 


D=D-[} a ok “P] FOR ALL 2 nt AND 0 <n, < ly 
1= 


ect, {n,@)} FROM D USING THEOREM 4.1 | 


f 


6] COMPUTE UNIT PULSE DIFFERENCE MATRIX: 
{h,,,(n)} = th(n)} - th,@)} FOR ALL O<nji < L 


eae, a 






AND MEAN SQUARE ERROR € 


OM 


10 
"| PRINTOUT __ 


H 


= 


FIG 9.1 EXTRACTION OF AN N-DIMENSIONAL 
ALL-POLE TRANSFER FUNCTION 


iS 





illustrate the application of the all-pole transfer function 
extraction method to a two-dimensional unit pulse response 
which 1s chosen for convenience to be normalized, i.e., 
h(0,0) = 1. The unit pulse difference matrix and the 
mean-squared-error are listed for every traversion of the 


iain loop. 


meee. step 1: (Notice that the step number corresponds 
to the block number in Fig 9.1). The given unit pulse 
response {h(n,,n2)} is listed in Table 9.1, where 
Mey = 10 , O < na < 10. The algorithm will be 


terminated if the mean squared error is less than ue” 2 : 


ai A graphical representation of f{h(n1i,n2)} 


mee. & = 10 
Mmeroewn im Fig 9.2. Step 2: Initialize D(zi,z2) , I , 
{a(x} Mecca Ali(n,,nod} = {h€my4mo)} .9 Steps 3,4:. The 


coefficients amen) and Te are extracted along the 


igen didgonal, i.e. n,tn, = 1 of {Ah,(n,,n2)} , where 


4(1.9) = +52 > and 


210) 1) eel 

The two-dimensional characteristic equation is identified 
mee D(z,.Z2) = 1- .52 Za - .65 25 - Step 5: The unit 
pulse response {h,(n,,;n,)} corresponding to 

MimeZ 1525) = Do +¢2,,25) is computed using the recursive 
technique of theorem 4.1. Step 6: The difference matrix 


{Ah,(n,,n,)} is obtained by matrix subtraction: 
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memes sn,)} = {h(n,,n,) - {h,(n,,n,)} and is listed in 


table 9.2. The mean square error ¢ 1s computed as follows: 
We 10 
1 a 
: 2a » yy Piva Gatien) =~ 2928) 
n,=0 n,=0 


7 


Meeps 7,8,9: The tests e« versus 10 ° and I versus 20 are 


traversed and the main loop index I incremented to 2. 


Meepee. oteps 3,4: The coefficients 


C20) = Nia 2 0.) = el 
acai el) = iat © es lh) a oa eZ 
(9.2) = 4h2(0,2) = 008 
Me@eerdentifiied off the second diagonal cut, i.e., n,; tn, = 2 ,j 
of {Ah.(n,,;n.)} . The two-dimensional characteristic 
equation becomes 
Mee.) = 1 ~ .52z, - .65z5 + .llz. + .22,] 23. - .003z3° 


Meeps 5,6: {h.(n,,;n,)} is calculated and {f{Ah,(n,,n.2)} 
Meee table 9.3), is formed. ¢ = 9.854E -3 . Step 7 to 9: 


Tests traversed and I incremented to I = 3 


Loop 3. Steps 3,4: The coefficients = Ah3(n,,n>2) 


a 
Ca 


for n,+t+n, = 3 are identified as 


Toe 





ais 0) = ate § 8) (0), = mony. 
2¢9.4) = Aha(2,1) = .0032 
aGue nen (0,3). = 0. 


The two-dimensional characteristic equation becomes 


D(z, .2,) = 1 ~ 0.522," - 0.652,' + .11zj° + .2z;°2)° - 0.0032," 
1 


ran | -~ = 
- 0.034z)° - 0.00322, 2) 


Beeps 7,8,0: Tests traversed, I = #4. 


Meee. «6 teps 3,4: The coefficients are read from the 


fourth diagonal of {Ah,(n,,;n,)} as 


2(4 50) = .027 


des } = -~.053 


mG — (22). = (054) 


The two-dimensicnal characteristic equation becomes 


ca | = ro 2 ee a8 
D(z,,Z2) = 1 - 0.522, - 0.652, + .1iz;’ + .2z; zp - 0.003z>- 


3 Bt ce} ail -1 «3 
wows e- 0M0S22, 2, = 0.0272, + 0.0532; zp 

Steps 5,6: The unit pulse difference matrix {Ahs(n,,n2)} 
Memeasted in table 9.5. e¢ = 5.07566E -3 


eps 7,8,0: Tests; I = 5. 
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Loop 5. Steps 3,4: The coefficients along diagonal five 


eee 

acy 7) 70-0351 

acy 3) = -0-0107 

WasnOne™ “Coe2) * “(15l) > 700.5) o 
and 


Meyz ) = 1- 0.522; - 0.6522) + .11z,° + .2z;'z2° - .00322° 
3 


= -~Z —1] —4 = ie 
—MOctzaw= © .00322; Zo = .027z2; + .053z) Zo 


~ 0.0351z) 22. + 0.01072) 25° 


steps 5,6: {fAh,(n,,n,.)} 1s listed in table 8.6. 


€ = [71-10 ¢ Sieps J,8,9: Tests: 1 -= 6. 
meoeeo: oteps 3,4: The coefficient Arg 9) = -0.067 is 
3 


the only non-zero entry along the sixth diagonal of 
{Ahg(n,,;n.)} . The two-dimensional characteristic equation 


becomes: 
a =i my. cr ae | =o) 
Meyez>) —- 1 - 0.522, - 0.65z2 + .llz, + .221 zo - 0.003z, 
3 a, ~4 i =3 
= 0.0342, = 0.00322, Lo = O02 7z24 3 OROS3SZ4 Zo 


SUGune sinus: + SGlldpe, zat 0.067z.- 


aS 


Steps 5,6: {Ah7(n,,n2)} is calculated. The mean squared 
error is identically zero. Step 7: The main loop is 


terminated by branching to Step 10. Step 10: PRINT OUT: 


ll 


H(Z1,Z2) = Se a es eT ey ssa 
imeoortemenn 652, + .llgge + .22, zo —.008z5° = .0c42;—" 


pa A PS 


EC ae r eee 02725 fF 05321 25. = .0e51z) z5 
ae —6 
aOLO IZ) Zo + .06/Z21 


Mete H(z1,Z22) is exact, i.e., the unit pulse response 
MieeH(21,22) is identical to the given unit pulse response 


{h o,n2)} 


D. EXTRACTION OF A RATIONAL TRANSFER FUNCTION 


An N-dimensional transfer function G(z) is defined 


» by zt 


= 2,(0< &.<L.) 
epee ee eee (9.2) 


in eq (A.10) as 


N 
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The Taylor expansion of G(z) is obtained by applying 
eq (4.21) in slightly modified form which takes the 


numerator palynomial into account: 


)) Sie. Co3) 


@(z) = 
n 
(O<n, <o) 
where 
g(n) = Y a- g(n-m) + B— (9.4) 
m © n : 
m 
Ooi, Se 
— 1-1 
N 
‘ Mm. ~% 0 
pel 


Equation (9.4) represents the recursion equation corresponding 
to eq (9.2) for a unit sample input (See Appendix A) and can 
be verified from Figs 3.12 and 3.13. The application of 


eq (9.4) is shown in the following short example. 


ee) 





Pxample 9.2: The Taylor expansion coefficients g(n, ,n2) 


for 


=i) el ail 
7 + Zz ar Zz + ZZ 
G(z)5Z>) = 


-} -) -)] -) -) -2 -2 -1 
Lt - 102. - Go1Z22 — 222, 22 - A122) 22° —- 2)2) Ze 


=i2 2 oa” (9.5) 


2041 ~~ Go9242 ~ Aog24) 42 


are computed with eq (9.4) to be 


g(0,0) = Boo 

(1,0) = a) 9g0,0) + By, 

meet = at BL 0,0) + By, 

ee = 8 24051) + ap, 261,90) + o,,800,0) + B,, : and 
Biels) = ,,8(n,-1,n.) + ap, g(n,.n,-l) + g,,g¢n,-1,n,-1) 


+ Op 22(N, ieee A) 7 By SIA =e 


Maenever n,; > 1, nz > 1 =~, and all initial conditions are 


zero. This result can easily be checked from Figs (3.8) 


mmmouch €3.11) in Chapter 3. 
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If fh(n)} represents the unit sample response we wish 
to synthesize by a recursive transfer function G(z) , as 
Metamed in eq (9.2), then by equating terms of {h(n)} and 


maen)} leads to the following equations involving the a's 


ema B's 
mon) = ) = g(n-m) + B= shor O<n. <b. a). [e) 
m 
Oem. <1. 
—ji- i 
N 
- m+ 21 (0) 
i=l 
and 
mr) = ) a= g(n-m) oor L.<n (9.7) 
m 
O<m.<M; 
meals? 
N 
» meer Oo 
a 
mi =( 


It 1s observed upon fixing L and M , defined as the 
total number of nonzero coefficients in the numerator and 


denominator, that there are not enough degrees of freedom 
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in G(z) to generate all {h(n)} , where 0 <n. < 


|; — 
sao 0 <i < N. &£ and M ean be written without loss 
of generality as 
N 
L = I CL; ame : C35 83) 
j=1 
and 
N 
Mee= I CM. te a (39.9) 
a=) 
For example, for eq (9.5) 
Bes = Clatl)Cbotl) + (Mitl)(Met1l) - 1 = 12 


which is obviously insufficient to generate an infinite by 
infinite matrix {h(ni,n2)} . Nevertheless, if §& is 
minimized, which is defined as the weighted squared difference 
between the given sample response {h(n)} and the 
approximating sample response {g(n)} up through the pth 


sample, where 


N 
no I (P; +1) Co e> 
1=1 


and P > L+M (9.11) 
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then 


_— a 2 
Ss )) WS ES - e(B) | ec 


where W— > 0 
p 


The minimization of &§ can only be achieved by iterative 
Means, Since g(n) is a nonlinear function of the a's and 


B's. There exists, however, one case: 
ot Ce esp) 


where the minimizing coefficients of eq (9.12) can be found 
by solving P linear equations! 


By assuming g(n) = h(n) where 0 < ns 


ae it is 
possible to solve eq (9.7) for the &'s that enforce 

Mem) = h(n) for L. <n; . With the knowledge of all a's, 
it 1s possible to solve eq (9.6) for the 8's, such that 


Me = h(n) for O <n. 


, £ Ll; . Equations (9.6) and (9.7) 


and the regions of summation are illustrated for the 


important two-dimensional case in the following example. 
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Pxample 9.3: Given {h(n,;,;n,)} , where 0 <n, < P, and 


mem, <P, , find 


Ly L, 
22> Z1 Ze 


R1=0 Qo=0 
SS > (9.13) 
= eg 100) Cities (CU, 
7 > ) Sh, 
mi=0 m2=0 
(m,+m2# 0) 


where p= [t+ MM. 
The denominator coefficients are computed from eq (9.7), 
where 


My Mo 


h(n,,n2) = > a Oram, ee a-M »sNa-Mo) Clea 


m,;=0O m,=0 
GitmsneD Ese 


The @'s are computed from eq (9.6), l.e., 


My Mo 


h(n, i>) = » » ne h(n,-m, »No-M,) t Ban, (9.15) 


m,=0 m,=0 
(m,+m,)#0 
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myeation (9.14%) is summed over all n, >L,, n,> L,, whereas 
mes. to) 1s summed over 0<n,slL,, and 0sn,sL,. This can 


be shown graphically in the following manner: 


Region of summation of 
eq (9.15) to compute a's . 
(L entries) 


Regicn of 


summation of eq(9.14) 


to compute 6's . (M entries) 





The above outlined procedure equates the truncated power 


series 
Gz) = » ING) wae (q16) 


to the first P terms of eq (A.11). It generates the first 
h(n) samples, where (0 < n; <P.) , exactly. 
The approximation of a power series by a rational 


mmetion 1S Known in the one-variable case as the Pade 


approximant. 


Example 9.4: Given a causal unit pulse response {h(n,,n2)} 
Bere 0 <n; < P} , O< ng < Pp , and B= C22), 
Extract G(z,,Z2) which is chosen such that the number of 


samples in the unit sample response equals the number of 


coefficients in the transfer function, l.e. 


-1 —1 


a 4 
+ + 
G(z,,z.) = Boo * Bio2 Bo1Z2 + B1121 22 es 


1 -) -1 -1 -) —2 - 
SeOrocwmeee Coveo 9 = 90174; °29 -— A292) ~ ApezZ2 


Comparing eqs (9.13) and (9.17) it is observed that 


Mee = Lo = l 
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Equation (9.14) becomes 


h(n, M2) = aioh(ny-1,n2) + aoih€ny,n2-1) + ay yh(ny-1,n2-1) 


at Qooh(n, yN2—-2) cf M2 9 (n4-2,n2) 


Momeni > ly or nz > lo , which can be written in matrix 


form as follows: 


Bie?’, 0 ) h(1,0) 0 0 0 h(0,0) 10 
2, 1) Biles 162.0) 7a 50) 0 nc Os) Oo 1 
BnG2 2) in@ie enc? I) nC) h(230) ht0,2) Cael (9.18) 
n(1,2) n(0,2) h(1,1) n(O,l) hQl,0) 0. Boe 
hed, 2) 0 Jal Oy), 0 enqeraep) 0 Q20 
Similarly, eq (9.15) becomes 
h(n, 5M) = ay 9h€ny-1,n2) + ag yh(n,5n2-1) + a, ,h(n;-1,n2-1) 
+ Qo ohn, 5M2—-2) + Q2 gh(n,-2 Mo) a ss 
Nyl> 
mor (0 SS 10) Ss 1 and OQ S No = HE 
or, 
h(0,0) = Bog (9.19a) 
mci,0) = aigh(0,0) + Bio (9.19b) 
mcO, 1) = &o91n(0,0) + Po CoE e) 
mic, 1) = a3 9n(0,1) 7 ao 1n(1,0) + a11n(0,0) + Bia (9 .19d) 
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The above computations are performed for a numerical 


example as follows: 


1 .053 mile 
mines) } =) .05 AS #005 
se .024 OR 


Equations (9.18) and (9.19) can be solved for a's and g's 


respectively. The result is: 


10 i... 0 SSO 
Oo1 cpl rel cae 
Qy = 5 EO 
9 2 Bplecnewls 
Q2 0 206559 
Bo 0 ag 

Bio -1.03807 
Bo. ~ 1.244 
Bif =) 201233 
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fee EXTRACTION OF A PARALLEL ARRANGEMENT OF LOW ORDER 
ALL-POLE OR RATIONAL TRANSFER FUNCTION 


It 1s important for the analysis, design and implementation 
of multi-dimensional recursive digital filters to approximate 
a given unit sample response or a high-order transfer function 
by a parallel arrangement of low order sections. It will be 
shown in the next paragraphs how both of the previously 
discussed extraction techniques can be applied in an 
algorithmic approach to extract a parallel arrangement of 
low order transfer functions. 

To approximate a given unit sample response by a parallel 
arrangement of low order sections, the following steps must 
be performed: 

- Partition a given unit sample response {h(n)]} , 


0 <n; 


[A 


I. in subsections each of size R , such that 


al 
mee ty 


7 


- Designate each subsection with the vector index j , 
imechie usual manner, where 0 < 5 <u; 
These two steps are illustrated in the following example 


for the important two-dimensional case: 


imeomple 9.5: The given unit sample response matrix 


{h(ni,n2)} where 0 <n, < 29 , and 0 < nz < 14 is 


Partitioned into R (4,9) size subsections, where J = (2,2) 


as shown in Fig 9.4 and designated as shown in Fig 9.5. 


It 1s noted that the unit sample response entries of the 


zth subsection are 


Rt 


hx(r) = {h€jiRa * ji + r1,---,3yR, + Jy + a) 
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for 0 See 


<< Re 
Ll — 


Gs CF eal.) CD Sine veoter  immae:: 


ietation, 
her) ee ing aeet 


The algorithm for the extraction of a parallel arrangement 
Sielow order transfer function, which is shown in Fig 9.6 
can be viewed using the above notation as an application of 
the all-pole transfer function extraction algorithm (see 
Section C), operating on the vector index Jj and on sets 
of unit pulse coefficients ha(P) 

It was shown in Chapter VIII that the unit sample response 
entry h(n) contains the corresponding transfer function 


coefficient aq) >» as well as non-linear terms composed 


of acy) » where ORD Sine) Sahn ei << any 
and Ip SS imp 
and ky < Ny > where 
N N 
oe 
al a 
= IL 1=1] 


(The notation used is the same as in theorem 4.1 and in 
Chapter VIII.) Also if the unit sample response of the 


Eransfer function H(z) having {arr} Soom ricienis . 
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Satisfying (i) < (n) is known, then a=) is equal to 


the qth 


entry of the difference between h(n) and the 
corresponding unit sample response entry of H(z) 
This concept was generalized in Section B of this 


Smapter to extract all Jan) , along the ma diagonal cut, 


where 
: : 
q hg = | ; 
ere 
at the same time. It is necessary at this point to further 


3 


generalize this method by extracting all sets of ay : ie eee 


ae 
a , where > ee ee a 
(r) if a — 


along the 1 oe 


diagonal cut in {he(>) } from the difference 
of f{h(n)} and the unit sample response corresponding to 
the sets of coefficients extracted along all F diagonal 
cuts, where F<1. 

This generalization follows directly from the propagation 


rules presented in the previous chapter, and will be illus- 


trated in the following example. 


eo: 


Example 9.6: Given a unit sample response {h(n,,n2)} 


of the size shown in Fig 9.4. 


Sa 32) Cu ad e) 
Beaxcract a and a 
(ry are) 0 Cane) aL 
Meee O ri = 5 , OS rp <9. The set of coefficients 
in 5 32) 
a | is extracted from f{h(n,,;n2)} using methods 
(rr; ,P2) 


discussed in Sections C and D, and the all-pole transfer 

fume ton age constructed. The corresponding unit sample 
response is computed and subtracted from {h(n,,n2)} 

The resulting {Aho(n;,n2)} has zero entries for 0 <n, < 4 


mons lo < 9 . The sets of coefficients 

om.) 2) (EAs, (G5) 
a = a ; a 

my , re ) (ri ,Y2) (ry ,P2) 

ak 

are read from {Aho (4tri,r2)} and {Aho(ri:,r2+10)} , 
for r;,Yr2 as above, using previously discussed extraction 
techniques. This procedure is written in algorithmic form 
for the N-dimensional case and is shown in Fig 9.6. 


The algorithm extracts a set of transfer functions 


{Hx(z) } which are arranged as follows: 


)) H>(Z) a) Conca 


J 
(for all j traversed 
in algorithm) 


H(z) 


ee 





correspond to a unit sample response which has exactly 
the same entries as the given unit sample response h(n) 


fem the first 


N 
I] SGRe Det 

a 1 1 
12H 


mmerres. ithe order of each transfer function H(z) is 


less than or equal to 


N 
> “i 
a 


aS Ih 


The following example is designed to illustrate the algorithm 


Mmaerag 9.6. 


meee a xXAMPLE 9.7 FOR E 

Given the unit sample response listed in table 9.7 and 
the termination criterion € = jo Es Gi hoes 
partitioned (step 1 in Fig 9.6) into subsections of 
Meeol.2) size and redesignated (Step 2) such that J is 
Memmeited as: uJ = (2,2) . It is noted that {h(n;,n2)} 
is not of square form and that the termination criterion in 


mock 10 (Fig 9.6) is therefore adjusted to I < 2 


mecordingly. 


weep 0: The main loop in Fig 9.6 is traversed for jitja=I1 , 


where I = 0. In step 2 (corresponds to block 2 of Fig 9.6), 


els 






GIVEN: UNIT SAMPLE RESPONSE {h(n) } 
ORDER R OF LOW ORDER SECTIONS 
TERMINATION CRITERION: MEAN SQUARE ERROR € 


oIPARTITION {h(n)} INTO R SIZE SECTIONS 
REDESIGNATE SECTIONS. IDENTIFY 7 








Saeed eR RE AR Pe eS Se ee 


3 
INITIALIZE INDEX OF DIAGONAL CUT I = 0 


ee : one. Se. 
TRAVERSE THROUGH ALL 3 3 ) 3, =I AND 0 <3; <J, 
a Ed 


; EXTRACT H(z) BY ALL-POLE OR | 


RATIONAL METHOD : 
asl al es 









eine COEFFICIENTS Hs ee 


SerPUrE UNIT SAMPLE RESPONSE hy Ci) Or 
H,(@) = J B22 GRD 


FOR ALL 4 TRAVERSED. 


ee ee re 





ee Oe re er Se ER Ee OO EE 


i 


: COMPUTE {Ah, (oe toca) | = thy (n)} 


AND MEAN SQUARE ERROR & 





lial 
{ it = singel E Ss 





FIG 9.6 EXTRACTION OF PARALLEL ARRANGEMENT OF LOW ORDER 


| ALL-POLE OR RATIONAL TRANSFER FUNCTION 
| een 


Hog 9) 62142) 1s extracted using the all-pole extraction 
technique of Section C. The resulting transfer function 
coefficients are listed in table 9.10a. In step 8 
Mene«n;,n2)} is computed, which is also listed in table 9.8. 


The mean square error 1s eé = .2226185 . The termination 


tests are traversed and I= 1]. 


Maem. is traversed for jitjo = 1 , Hey) 621222? and 

Hg 1) 621922) are extracted. The coefficients are listed 

in table 9.10b,c. The unit sample difference matrix of 

step 8 is shown in table 9.8. e = .1558379. Tests traversed 


ama tf = 2 


Meep 2: jitjo = 2 the coefficients of Hog 9) 621922) 

Ho, 1) 621 222? and Wg a SMa eee are extracted and listed 
in table 9.10d,e,f. The unit sample difference matrix 
{Ah2.(n,;,;n2)} = {0} and mean square error test (step 9) fails, 
S@aethe algorithm is terminated. The printout in step 12 of 
the extracted parallel arrangement of third-order transfer 


Mmemetron {H,. . .(2;,;22)} is as follows: 
Cee) 


-2 

mee Z 422) = Hg 9) 621222) + Hg) 621922921 
ap dat Ge za elas C25 eee Pees 

COs 2 (eibesis) 2 


= =o 
+ Heo 9) 641222) 21 + Hog 2) 621222) 22 


Note that the test in step 8 will only be utilized as an 


additional regulative to terminate the main loop traversion. 


21S 


G. DISCUSSION 

There were essentially three methods developed in the 
preceding sections to extract an all-pole, rational and 
parallel arrangement of low order sections. The order of 


an all-pole transfer function extracted from f{h(n)} , 


0 <0 <I; temas Loulows 
N 
order H(z) aS q I. 
= |) 


The rational transfer function method extracts a transfer 
function of order determined by the choice of L and M , 
mmere £ + M equals the number of entries in {h(n)} 

The low order all-pole or rational extraction technique is 
an important result by itself, since it represents an 
approximate partial fraction expansion technique in 
N-variables which cannot be achieved exactly. (See 


montactorability proof in Appendix B.) 
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Poot GN OF N-DIMENSTONAL RECURSHVE FILTER 


io «60OLNTRODUCTION 

The time domain design of recursive digital filters is 
achieved by selection of transfer function coefficients such 
that the corresponding unit sample response approximates 
a given one and by checking the stability of the resulting 
transfer function. Both steps of the time domain design 
procedure have been developed: the selection of coefficients 
by extraction techniques in Chapters VIII and IX, methods to 
analyze stability in Chapters V, VI, and VII. At this point 
it 1s necessary to establish a measure of quality for the 
approximation step. All approximation procedures produce a 
[meeemer function HCZ), such that the associated unit pulse 
respone has identically the same entries as the given unit 


sample response {h(n)} , where 0 <n. sore the 


ee 
a 


1 
first k terms, where k; < I. . The quality of approxima- 
tion will be identified by comparing the magnitude squared 
spectra of the given unit sample response to the magnitude 
squared spectra of the extracted transfer function. 

In the following section three design examples are given 
omirdstrate this measure of quality. In Section C several 
aspects of frequency domain design are discussed and the 


applicability of stability analysis as well as approximate 


partial fraction expansion techniques are pointed out. 


230 


Finally, areas for further research based on the presented 
techniques are outlined in the summary at the end of this 


eiapter. 


fee tT MeE-DOMAIN DESIGN EXAMPLES 
Peeoxample 10.1: All=Pole Transfer Function 
From a given unit sample response {f{h(n)} , 


0 =) 1. 


+ <I. » an all-pole transfer function is extracted 


using the algorithm of Section C in Chapter IX (Fig 9.1). 
The order of the approximating H(z) is chosen by comparing 
the quality of approximation in the frequency domain. The 
stability of the resulting recursive digital filter is 
analyzed with conditions given in Chapters V to VII. 

Given a unit sample response {h(n,,n2)} , 
C0 Sy Rezo. 5) 2h 15 observed EMat the magnitude squared 
spectrum is symmetric, i.e. JHC, 50) |" = JHCO, | ° for 
meieew , Which is shown in Fig 10.1 Cor 0 eu Se) s 


p 
The application of the all-pole extraction technique for 


i < ns Peeetcddaeto a Didimear transfer function (21.22), 
where 
E: dk 
H(z21,22) = ——————_———————— : Cio ay 
ay -1 -1 -1 
aro eee Coi2o ) «e112, Zo 


Meemiransfer function coefficients are listed in table 10.la. 
tjijgenextraction of a fourth-order transfer function of the 


form 


2S 





i 


L = 1021 - 0122 -— 112) Z2 — A202) 


Seon a1e Ze = OooZt = O22 22 


ee: aco Co a ClO sy) 
SeOleeilaeo = O99Z2) S22 


(Coefficients in Table 10.1b.) The corresponding magnitude 
squared spectra are observed to be symmetric and the 
comparison in Fig 10.1 along the w,-axis shows that the 
fourth-order approximation has good quality. The transfer 
Met 1Ons in eq (10.1), (10.2) are identified to be stable 
by Tables 7 and 6.1, respectively. 

MeeexomplemuOn2; Rational Transfer Funetion 


Bommemeavensencn) 95 O <n. <0, 5 @ rational 


transfer function is formed by solving 


P = | (P.+1) 


Bigecam equations, where 0 < P. ama f = aah 


ele 
ieee 


al 3 
as defined in eqs (9.8), (9.9). 
For a given th(ni,n2)} , O<n; < 25 the associated 


symmetric magnitude squared spectrum is shown in Fig 10.2. 


The second order approximation 


(Sz 


Boo + Sana pneu uy Oo eA (10.3) 


a a | 


FOZ 1 522) = a ae < 
Pas unas = ey ao = Oy yey iz 


and the fourth order approximation by 


me = oe -2 
Scommeciez: ' bots  tebrizi co. + Booz) 


H(Z14Z2) = = -1 -1.-1 = 
= @192: -— 122 — G112) 22 = Go02Z1 


-—2 -2_~1 -1 -2 —2_=2 
teogiezo. + ome Zo Hii Z) Zo + Boo21 Ze (10.4) 
= are ae ee 0 aD : 
meee 22121 -242 = Alot Zo = G2221. Zo 


for which the coefficients are listed in table 10.2, result 
in symmetric spectra, which are shown in Fig. 10.2. Again 
the fourth order transfer function gives a good approximation 
The stability of the two transfer functions is identified 

by applying Tables 7 and 6.1, respectively. 


o. Example 10.3: Parallel Arrangement of 
Low Order Transfer Functions 


The transfer function extraction is achieved using 


'the algorithms in Fig 9.6 as well as the all-pole extraction 


algorithm or rational extraction method described previously. 
The given unit sample response is partitioned and redesignated, 
as explained in section E, Chapter IX. For the unit sample 

response given in example 10.1 the partitioning is made into 


bilinear sections. The second order approximation is 


ae 8 


necessarily identical to the one in example 10.1. The 


approximation by three bilinear sections, i.e., 


— 2 —2 


for which the coefficients are listed in Table 10.3a, has 
a magnitude squared spectrum shown in Fig 10.3. The 
stability of each section is verified by application of 
Table 7. 

It 1s noted that the approximation by one fourth 
order all-pole section, as derived in example 10.1, exhibits 
a much better quality than the approximation by three 
second-order parallel transfer functions. 

4, Discussion 

All of the discussed design techniques do not lead 
to a locally optimal solution, as an iterative technique 
will, but they produce a viable solution in a comparatively 
short time. It is noted that the rational transfer function 
extraction technique for a numerator polynomial equal to 
a constant and the all-pole transfer function extraction 
algorithm produce identical results. 

C. COMMENTS ON FREQUENCY-DOMAIN DESIGN TECHNIQUES 

EOk RECURSIVE DIGITAL FILTER 

The availability of a general method to derive stability 
conditions in N-dimensions will remove one of the major 


obstacles to solve frequency domain design problems. Also, 
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the possibility of approximating a high order transfer 
function by a parallel arrangement of low order functions 
allows the realization of complicated systems similar to 


mre use Of canonical Sections in one dimension. 


D. SUMMARY 

The N-dimensional time-domain design techniques represent 
a step towards the use of recursive techniques to perform 
multi-dimensional filtering. Based on the theory of 
structures, stability analysis and design methods, future 
research is suggested to develop frequency domain design 
techniques, to generalize Parker, Girard, Souchon's [42] 
study of effects of correlation and structure on quantization 
Meese tO N-dimensions, and to formulate conditions for the 
existence and bounds of limit cycles in N-dimensional 


recursive digital filters. 
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fae ejay IES a\ 


See ae Oe DEG TTA sFIiLTER THEORY 


N-dimensional linear digital filtering can be 


described in terms of N-dimensional 7-transforms, 


Pe DEFINITIONS 
1. N-Dimensional Signal 
Consider the signal x. to be a function of N 
mmerab les, 1.é., TisseeoTys which span a N-dimensional 


space. <A uniform mesh on this space is chosen 


such that the coordinates of any Vertex ofthe mesh will be 
(niTi,n2T2,...,ny7T,) , where Tr = (T1,...,T,) are the mesh 
Spacing along the Tiye+eTy axes, respectively. Any vertex 


can be identified by the coordinates 


cn) CM ees hee ; CAR 
ei@ere (1) = (1) and the signal at that knot is denoted 
as x(n). {x(n)} is defined to be an N-dimensional 
mecauence, where each -© < n; Ke 
2. Linear, Shift-Invariant, Causal System 
A system is defined mathematically by a unique 


transformation or operator P that maps an input sequence 


fom) into an output sequence y(n) , i.e., 


WAGhiae= xno] 
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An N-dimensional system is linear if the principle 


of superposition holds, i.e., 


BileearGm)) + bx>(m)} = aP[x,(m)] + bP[x>(m)] 


A shift-invariant system is characterized by the 
Bioperty that if yn) is the response of x(n) , then 


y(n-m) is the response of x(n-m) , where 


(n—-m) = (n,-m,, no-m., ..., Ny My? 


A system is realizable or causal if changes in the 
output do not precede changes in the input. If Uo (n) is 


the N-dimensional unit sample defined by 


1 if eS ST Sa i ee 
Dion) = CAD) 


0 otherwise, 


then the unit sample response, 1.ée., 


Gnome =P CUgGn)) 


of a linear, shift-invariant causal system is zero whenever 


at least one component of (n) is less than zero. 
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3. SComvolution 
A linear shift-invariant system* is completely 
characterized by its unit sample response h(n) , ice., 


the response to the input x(m) is 


om) = » x(m)h(n-m) Y h(m)x(n-m) (A.3) 


m m 
~o<m.<o —©<JIL. <0 
C3 aS ) ( <M. < ) 


where 


it! 


ee a ae 


Equation (A.3) defines N-dimensional convolution and will 
Bliso be represented by y(n) = x(n) « h(n) 

The difficulty of applying eq (A.3) in two-dimensions 
in all but the simplest cases is illustrated by Rabiner and 
Gold [2] and can be interpreted as an omen of some of the 


inherent problems in working with multidimensional systems. 


%e 


For causal system the lower limit is zero. 
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A system is said to be separable if its unit pulse 
response h(n) can be factored into a product of one 


dimensional sequences, 1.ée., 


h(n) 


it 
— 
Di 
bu. 
oN 
= 
te 
/ 


(A.5a) 


The advantage of a separable filter is that the multidimensional 
convolution for a separable input can be carried out as a 
sequence of one dimensional convolutions. This follows from 


rewriting eq (A.3), for 


N 
incm)) = II h; Gn; ) 
Is I 
and 
N 
cn. = II x, (n;) 
WS ih 


245 


as 


N N 
vin) = » I h, (m; ) II x. (n;-m;) 
ene. Nera ew 
which equals 
= > hi€m,)x1(ni-m)) yy heo(m2)x2(n2-m2) 
m)=--@° Mo x==—-O 


y hy (my) xqy(ny-m) | qve.d. — (A.SD) 


Bs Z-TRANS FORM 
fee betinition 
ine N=dimensionaL two=Sided Z-transform X(z) , 


mere (7) = (21,22,.. ) , of the N-dimensional input 


Zn 


sequence x(n) is defined by 


Mem ee eZ(x(n)) = » Gnas iz (A.6) 


n 
=—o<n. <0 
(-<n, <e) 
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— =p O05) = Ie> =i) 
0) N ; 
where Z = 4F ez} a alge ilove Z- S are complex 


variables, which, expressed in polar form, 


th 


where Ws 1S timer 1 Spatial frequency variable and can be 


used to rewrite eq (A.6) as 


X(Z) = »» (A) ( ee ot (Ne 78) 


nN 


Similarly to the one-dimensional case, the N-dimensicnal 
Z-transform can be interpreted as the N-dimensional Fourier 


transform of x(n) , whenever 


For convergence of the N-dimensional Z-transform, the 


sequence x(n) z must be absolutely summable. This is 


equivalent to 


)) |x(n)r-" | < @ (A.7) 


ro 
(-o<n. <o) 
a ee 


2U7 


The set of z-. , for which eq (A.7) holds, defines the 


iL 
region of convergence. Notice that because of the multi- 
m@ication of the sequence x(n) by rv , it is possible 


for the Z-transform to converge, even if the Fourier 
transform does not. 

Cee CODerties 

Properties of the N-dimensional Z-transform which 

are needed in subsequent chapters are derived as extensions 
of the well-known one- and two-dimensional cases. The 
proofs are omitted whenever they follow directly from the 
definition of the N-dimensional Z-transform. The properties 
of the N-dimensional Z-transform are summarized in Table A 


a. Linearity 


miieex(m) + b y(m)] = a Z[x(m)] + b Z[y(m)] 


where a and b are constants. 
IS Sigal seme ye 


If ym) = x(m-p) 


then 


Wz) =z? Xz) 


ee Cem 6 Ubon 
If y(n) is equal to the convolution of two 
N-dimensional sequences x(n) and h(n) , then the 
Z-transform of y(n) is equal to the product of the 


@-transforms of x(n) and h(n) , i.e., 
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ci 


) x(m) h(n-m) 


m 


then 


WCE) = Co Im) 


This can easily be seen by writing 


LO) Y | Y x(m) h(n-m) 3S 
n m 


fresny so) (ccm; co) 


and after interchanging the order of summation and changing 


the indices of Summation, this becomes 


we) = y Ges Y n(k) 27k { i-m 


k 


(-e<in, <co ) (-e<k, <oo) 


ea) 
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which is equal to 
aay =) GZ Cz) 


3. imverse Z—-lransrorm 


MictiiversiOm LOmmullastor the vabeweademens 1! Omnia s 


given by 
Gm). = a pe. f X(z) zsmet) dz (A.8) 


— 


where the N paths of integration Cy are within the regions 


of convergence of eq (A.8) and dz = dzi*dz2‘**dzy This 


fermula tollows by substituting X(€z) in eq (A.6) , which 


becomes 


ac) = + Meg db x(n) 20 pera ae 
(213) J f 


Q 
| 


The interchange of summations and integrations, which is 
justified by the absolute convergence of the series for 


7), yields 
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By the Cauchy Integral Theorem of N-variables and since 
the paths of integration Cy are chosen within the region 


of convergence, it follows that 


(any) if Cen) = (5) 


— 0 otherwise 


Pee RECURSIVE EQUATION 

An important subclass of linear, shift-invariant 
N-dimensional filters is constituted of those systems for 
which the input and output satisfy a linear, constant- 


coefficient difference equation of the form: 


ZS At 


N 
where ) ie ae) 
i 


with aca 


and Gives © Cla. . 0 


are the sets of constant coefficients which characterize 
a particular filter. 

Equation (A.9) does not uniquely specify the input-output 
relationship of a system. This is a consequence of the fact 
that as with differential equations, a family of solutions 
may exist. Therefore, a set of initial conditions must also 
be specified. 

It 1s assumed that if a system satisfies a linear- 
constant coefficient N-dimensional difference equation, it 
will also be shift-invariant. 


Pmemene set of all a(n) 1s non-empty for all m; such 


. Mere 0 5 
a 


1=1 


telat 


then eq (A.9) characterizes a recursive system algorithm and 
may be used as a computational realization of the system, 
either by programming a general purpose digital computer or 


by implementation using special purpose hardware. 


Boe 


D. TRANSFER FUNCTION 


Taking the Z-transform of the N-dimensional linear 


difference equation, eq (A.9), leads to 


NCA) (Co mm (ee 





where 
2% 
») DE) z 
gv 
= Chit <ias) 
H(z) = =a (A.10) 
, -Mm 
il - a an) a 
ay 
‘aa 7% 0) 
v—s)) 
N 
Note » m. = OP aneinvess T= 0), > Mm, = 0 
See 


Equation (A.10) is shown in Table A.2 in long hand notation. 


The order of H(z) is defined as Pree unien sald M. : 


me. , 


N 


Omder cz) = yy M. (A.lla) 


1=1 
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The degree of H(z) is defined as follows: 


Degmee H(Z)"= max (M. ) CAS Iai) 
i 
all 1 


A coefficient Ars) of the N-variable characteristic equation 
1s called a leading coefficient, if its order equals the 
order of H(z) . Note that there may be several leading 
coefficients. 

For H(z) as defined in eq (A.10) with constant real 
Coefficient and for a causal unit response sequence h(n) , 


with real entries, 


0<2.<L, ie 
2S Sn Y 116s) ae Ae (A.12) 


where the right hand side is the Taylor expansion of H(z). 
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Pee OLEARABLE TRANSFER FUNCTION 
For a separable system, as defined in Section A.4, the 
transfer function H(z) ean be written for a separable 


input as 


This is proven by taking the Z-transform of eq (A.5b). 


ENS: 
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Ue le lagi) DIE. 9) 


PROOF THAT MULTI-VARIABLE POLYNOMIALS 
CANNGT BE FACTOR ED IN GENERAL 


To establish that multi-variable polynomials cannot be 
factored in general, it is sufficient to show the insolva- 
bility of the problem in one particular case [54]. If we 
assume z,° + z,° + 1 = 0 can be solved into a product of fac- 


Beiesewhich are integral in z, and z,, 1.e., 


1 


(pz, +qz, +r) (p' z,+q'z, +r" ) 


N 
Nn 
an 
N 
) 
ae 
= 
HW 


Rc MnCl Cua aan, 
* (pq'+ p'q)z,z, 

TP QO aca jb) Je!) van 
(Capa es ) (169) 745 


Since this is by hypothesis an identity, we can equate coef- 


mrelents and obtain: 


pp' = l (B.1) pm a jon @ = 6 (B.4) 
gigi! =") 1 Cee) jeuak ae ig" 0) CBE os) 
mre = | GB.3) gio Gl se = [ (B.6) 


First we observe that, on account of equations (B.1), (8.2), 


mums. 3), mone of the six quantities, p, q, r, p', q', r' 


gem be zero; and also p’' = = q' = Series Winslsi, cls el 


2) ae Nn 
logical consequence of our hypothesis, we have, from (B.4). 


5-5) and (B.6): 
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P Gd ig 

pee ees Bad 
q + 5 ( ) 
re ea (B.8) 
Goa 

Fee = 0 GE) 
r+q 


and, from these again, if we multiply by pq, rp, and qr, 


respectively, we get 


aan G- =a (B.10) 
D> ay 1 eaualG | (B.11) 
Geese) = 0 Cee) 


Diceerca@ = 0 (B.13) 


And from addition of egs(B.10) and (B.13), 


Po) = 0 
Beotm@mrents it follows that p = 0, which is in contradiction 
with eq (B.1). Therefore, the solution in this case is 


Mmpessible. q.e.d. 


2S) 


Jategle a) IDE). 16 
Peet OF GONDITPONS FOR BIBO STABILITY: 


SE = )» In(n)| < Cem 


N 
-O<N.<@ 
(-<n, <@) 


Prove that, 


1s a necessary and sufficient condition for BIBO stability. 
ieeeq (C-i) is true and if for a bounded input x(n), l1-e., 


een) | <M for all (n), then 


acm |=] Y) h(m) x(fi-m)| < 


Mm 
—-o<cm. <@ 
( oe 


<M )» lh(m)| < (Cae) 


mn 
(-o<m. <o) 
oe 


Thus the output sequence is bounded. The converse is proved 
by showing that if S = ~, then a bounded input can be found 
that will cause an unbounded output. Such an input is the 


sequence: 


1 af h(-m) > 0 
on sa= 


-l1 if h(-m) <0 


*This extends a proof in [2] to N-dimensions. 
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Baomeeq (A23), the output at (n) = (0) is 


WACO)! 7= ) (Gm) tole) = ’» [h(-m) | 


m mM, 
(-o<¢m. <o) (-o<m. <~) 
5 iiss ca ics 


@Wh2eh is equal to 


y(0) = » In(m)| = 8 


m 
(-w<m. <o) 
ae alg 


Fherefore, if S = », it follows that the output sequence is 


unbounded. 
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APPENDIX D 


Vv 
ULTRASPHERICAL POLYNOMIAL C0) 


The following paragraphs contain the derivation of some 
important characteristics of the Ultraspherical polynomials, 
which are based on several papers published by L. Gegenbauer 
before the turn of the century. 

Gegenbauer defines the Ultraspherical polynomial CMC x) 


in the following manner [62]: 


v ain (ntv-1)! 

ee 2 TT ant 
x 1 0 0 0 0 

1 
a(VFL) x gE ¢ 0 0 
Zita 
7 e 28 @ 
ae OGG 8 1 { 
Oped) 

0 0 0 0 me a: 
0 0 0 0 (n-=1) (n+2v-2) . 


UCnty-2) (ntv-1 
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A. Even/Odd Function Characteristics 
Gegenbauer develops in [63] the following relationship 


of Ultraspherical polynomials: 
v _ fa v a v 2 
be CX) 2(nty-1) x C__ jx) + Cnt2y-2) Cl_4(x) 0 Cie) 


Differentiating OnE) with respect to x , by differen- 
tiating eq (D.1) and by writing the resulting sum of the 


n determinants, as a determinant of Gene order leads to 


vtl 


Sela oe CDs 


d Vv 7 
A Ic (x) J = oe 


Similarly the pth derivative becomes 


dv 
gle 


r (rtv-1)! VIE Oy) 


Vv = 
[Co (x) J 3 9 sSDE ahi 


If we replace x by (-x) in eq (D.1), then the determinant 
is multiplied by ea, Therefore the Uitraspherical 
polynomial exhibits even/odd function behavior for n 


even/odd respectively. 


Pemberivation of Second Order Differential Equation 


to which oC) 1S a Partwcular ssoluctaon 


Gegenbauer lists the following values of Ultraspherical 
polynomials evaluated at x = 0,-1,1 , which are needed for 


the following derivation: 
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Vv = 
rma 0}) = 0 


e7C1) = Cnet 
n (2y-1)! n! 


(D.4) 


Qe) = C= 1) sep 
nh 1@) 


n Gi 
(v-1)! n! 


0) C=1) 


2n 


Gegenbauer shows in [63] that the following set of 
formulas can be found: 
C2 


vue Wor db = Wea 
Bex) — © _54(x) == Tee 


Vagal ate 


pee x) - CO (x) = = C(x) 


which can be written as 


1G Cx] 


ni 2G) tema) 


rev _ a 
(Ca een an 


1 a 


[c¥(x)] - Stc¥ (x) 


Vv 
ie n Cox) 


ae 
m 


and are used to prove the following equations: 


Walk 


2 
mt 1—-x* ) Ch-1 


(x) +/nx CY¥(x) - (nt2v-1) Co_j(x) = 0 (D.5 ) 
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and 


uv(vtl)(1-x*) Cae.) - 2v02vtl1)x ee) 


+ n(nt2v) C(x) = (D.6) 


From eq (D.6) 18 is observed that the Ultraspherical 
polynomial is a particular solution to the second order, 


linear differential equation: 
(l-x*)y" - (2vtl)xy' + n(nt2v)y = 0 OB) 


The identification of the differential equation (D.7) allows 
mars to relate Cc. (x) to other polynomials. By choosing 
v=0 eq (D.7) 18 solved by Chebyshev polynomials of the 
second kind. It is also observed that eq (D.7) is a 


degenerate case, 1.e., a: $8 of the differential equation: 
(l-x*)y" + [B -a - CatBt+2)x]y' + n(ntatBtl)y = 0 


which 1s solved by Jacobi polynomials. 


ere Roots of OMG) 


It can be shown that the roots of C(x) are real 
and distinct by assuming x= x; is a multiple root of the 


ultraspherical polynomial, which implies that the derivative 


SV arall 
n-1 


pms tO zero if and only if Com, 1s zero. Following the 


Bec (x) , i.e., C _5(x) must also vanish. Equation (D.5) 
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same line of reasoning it can be shown that x= x, must 
Pomoc or dll Ultraspherical polynomials, i-e.,; of 
Bee) Co 4x), .-., Cox) . But since Cj(x) = 1 it is 
concluded that Ultraspherical polynomials have no common 
roots and alo) has no multiple roots. It is easily 
identified from eq (D.4) that x= +1 are not roots of OC). 
ie sbunetLons = Cy <x) and Cro) have equal 


Signs at the roots of Cy (x) provided x is restricted to 


the interval -1l to 1; [65] . Gegenbauer treats the 
monet Lon 

Vv v v 

Ci Ox), C71) > Marrs Cin aR) 


as a Sturm sequence which has, by eq (D.4), at x= 1 the 
same Signs, and at x = -l1 alternating signs. 

It is therefore concluded that Ch (x) has n_ real 
pi@marstancet roots in the interval -l, ..., tl 


D. Representation of CACO as a Finite Sum 


Equation (D.2) combined with the fact that CGo ue el 


can be shown [66] to lead to 


19 G5 9 Vaust?2a © eae a 


v : 
Ciftl) = C21) Sp—eagrayr 


which is a special case of 
fim 2] 


eo ) S = (si) (ntv-A-l1)! (2%) 2724 
n\* mee a. Ee TSN GS IDs 


A=0 


In Fig D.1 and Fig D.2 the Ultraspherical polynomial 
ex) tome@aphed fom 1) = 58 and wv = .2022)1 “and 


n= 2(1)5 and v= 5 respectively. 
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APPENDIX E 


MEMORY REQUIREMENTS FOR TWO-DIMENSIONAL 
REGURoi VE, VRGl tA: Tilak 


Mitra et al. [67] have shown that the storage requirement 
of a two-dimensional digital filter realization does not only 
depend on the number of delays but also on the input output 
array sizes and on the computation sequence. This is best 
illustrated in the following example: 

Consider the two-dimensional fourth order coupled filter 
block diagram in direct form 3 realization, as shown in 
feo. 14, for Gde59 = B20 = 0. Figure E shows the same 
realization in a way which graphically illustrates the 
computation and storage required. If the output is calcu- 
lated row by row, the input must be read in row by row 
and elements of the preceding row of w's (defined in eq (1.2)) 
must be stored in order to calculate each new w value and 
hence each new output. After the w values in the preceding 
row have been used for the last time, they may be replaced 
by stored values of the new row. Therefore, memory required 
Will be Mt+t2 elements. If the output is computed column by 


column, then the memory required is 2Nt+tl elements. 
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APPENDIX F 


Hono 10 DETERMINE CANONIC STRUCTURES POR 
FULL QUADRATIC TRANSFER FUNCTION IN TWO DIMENSIONS 


Parker and Hess published in [38] a procedure to derive 
canonic sections in one-dimension. This method can be 
generalized to the multidimensional case. The immense task 
of deriving canonic sections in N-dimensions is beyond the 
scope of this study. However, a method to derive canonic 


sections for a given full quadratic two-dimensional case, 


mye. , 
| = =) =e) =? =2 
emis) | ¢ Poize + BayZp Zo ot Booze) booze 
H(z 5Z2) = <<, Seo) aan | oe Gaye) 
ieewctro2)  ~— C6122 = 03321 Zo = "@en4i = Oooz> 
is outiined in the following paragraphs. It consists of a 


two-dimensional state equation formulation using tensor 
notation for the full quadratic case, which is then 

Z-~transformed and solved for the underlying transfer function. 
A coefficient comparison with the transfer function (F.1) 


can then be used to develop canonic sections. 


A. STATE EQUATION FORMULATION 
The general set of state equations describing the one- 
dimensional digital filter, in which the state variables are 


labeled X , the subscript identifying the state variable 


- 


Zife 


and the argument its Value ata Particular instant Ot 


discrete time, has the form {31]: 


Ga) = A X(n~1) PCr) (Ph S2) 
Where A is a square coefficient matrix and 
_ a 
X(n) = [x1 (n) 2) een |] (F. 3a) 
and 
DO SENG nem) 84, ye (F. 3b) 
are the state and input vectors, respectively. 
If we extend the notation of equation (FG) mto eo 
dimensions we obtain: 
X(n} nz) = ((A)) X(ny-1,n) + ((B)) X(n; ,n2g~1) ig U(n, Mo) Cras 
where 
X11(n) nz) X12(n} 4n2) X13(n} 4n2) 
X21(n} ,n2) X22(n} 4n2) 
X(nj sno) = 
ayo) ° 
and 


edt S 





Uy, 1, (nq sN2) Uy 2 (nq »N2) U, 3(Nq »N2) 


Mo (nysNs? Uso iysnloy 
WG, sn) = 
TERE GohET ere, 
emaee<((A)) , (€((B)) are coefficient tensors. For the two- 


dimensional second order case, ((A)) can be written in 


ieensor preduct form [33] as: 


((A)) = 


[=I 
x 
[=I 


where M and N are the component matrices of Tensor A. 


But 
Mi, Mio 
N,, Noo No, Noo 
M x N = 
Nyy Ni2 Ni Nio 
Mo, Moo 
No, N22 No. Noo 


Pemilarly, ((B)) = Pox @ : 
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If the procedure of matrix elongation [68] 1s applied 


to the state variable and input matrices, eq (F.4) can be 


rewritten for the second order case in the following matrix 


form: 
Xy 1 (nq yng) Nii Ny» 
Mi Mio 
X12 (nq M2) No, No 


Bee KN) sT5 ) 


Weal, m7 | 


Qir Qie 
ai Pio 
Qo1 Q 
Qi1 Q)2 11 
Po, Poo 
21 Qs 21 


X12 (n, >No-1) 


| 


had 


X11(n, sNo-1) 


X09 (n, sMo-1) 


Xy 1(ny-1,n2) 
X12(n,-1,n2) 
- X,,(n,-1,n,) 
No | NG) Gs allpoe yy 
U,,(n, sno) 


Wi On Np) 
Gee) 


Xo,(n, sMo-1) Up ,(n, Fal) 


Us 2 (ny glo ) | 


The output equation corresponding to (F.3) has, for the 


case of single output, the form 


X,,¢n,-1,n,) 


X,o(n,-1,n,) 


y(n, ,n2) =e Ro R3 R, | 


X, ,(n,-1,n,) 


Xoo(n,-l1,n,) 


U, ,{n, Sil?) 


Ula gala eae ) 


Uz 1 (nq Nz) 


|U22(n, 5M) 


2D 


+ [S; S, 


Nai Ce ele) 


X12 (n, »n,-1) 


Si] 


D3 
X, ,(n, »n.-L) 


Xoo (n, sNg-l) 


(F.6) 


ine Z-transform of eq (F.5) and eq (F.6) is solved for 


NG ARR AOD, 


. i aeser®) 


Peeeetficient comparison of H(z1,22) as defined ain eq (F.1) 
and H (21522) leads to a set of equations with many 
solutions. However, if the Parker-Hess rules [38] generalized 
to two-dimensions are applied, then the above equations can 

be solved for a unique set of coefficient tensors ((A)) and 
((B)), which in turn can be used to construct a set of 


unique canonic sections. 
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APPENDIX G 


LOW ORDER TWO-DIMENSIONAL STRUCTURES 


Direct form and canonical structures of the third and 
fourth order uncoupled, as well as second and fourth order 
coupled, transfer functions are presented for the important 


two-dimensional case. 
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APPENDIX i 


COORDINATE TRANSFORMATION IN TWO-DIMENSIONS 


‘ RECTANGULAR TO DIAGONAL COORDINATE SYSTEM 


Given {h(n,,;n,)} ,; which has the form 





eGo sn (0) MCOnl) GO») 


GUO pbemCie ly) Inne @> 


EG) = nC), 0) Ln on oe) 


The transformation to the diagonal coordinate system 


th ¢m, >m,) } 5 where 


elo) In (252) 





ERO) yn 82) 


ievemy,.im>)} = 
gy We guile aA) WaGoil) 62) 


is achieved by setting 


mnee= Mi tno); and 


Mo = n2 


CoS 


Example: 


from 


EO 


H (2122) 


Be 
where 23 


=) 
and 22 


ime bilinear trancrer tUnctioen to etranctormmed 


2 A 
= = apa 
UF O59 2) ee Colo eZ ee 
2 IL 
Ss eee une ee | 
lL - 413031 —- 0181 22 —- 411321 82 


delays exactly as 21 does, 


delays along the diagonal entries. 


294 





APPENDIX 1 
Piel lel MaMUO NSC a Sle La seelvis Iu stOD) 


rOhetne DERE ALLO! Or 
STABILITY CONDITIONS IN TWO-DIMENSIONS 


Pela kRD ORDER CUNCOUPLED) FILTER 


im Unt Sample Respense 


The planar unit sample response A 


will be computed columnwise using eq (4.21). 


meme n>, = 0 (Column 0) 
imeoe Oo) = 1 
pe, 0) = G10 
2 
me2.0) = %10 + QoQ 
3 
MS) = aio f 20,0820 
k-2 k-2)(k-3 Kk 2 
ht(k,0) = a amc 1 Og © Cede eS 2D A10%29 
k-3) (k-4)(k-5) k-6 3 
J ( Dae a Aig 20 t+ eer + Eg 610 2%20) 
k/2 
A299 
a ee (k-1)/2 
where Eg § G10 2A20) = Oo) Oy -1 
ie I %10%20 
ELLEN Belge Tr PEI (oe ( 


aD 





The pn term h(€k,0) can be shown to be equal to 
‘{k/2} 

mdc,0) = » om heme) (k-2mt+2) pe (I. 
m=0 

which is equal to 
[k/2] 

Mik.0) = Geom oe a Cig2 
m=0 

Bemeene = 1 (Column 2) 

h(0,1) = ao, 


hHC1,1) = 2091010 
iG?) = ccen + 2091026 
3 
mes, 1) = Yap 1010 = 6a 1019020 
4 2 2 
HC4,1) = Sag 0.0 + L2apidioa2o + 3091020 


5 2 2 
hC5,1) = 6001010 + 20001010029 + 12091010020 


k(k-1) . k-2 if (k-1) (k-2) (k-3) _k-4 2 


h(k,1) = (kt) ay sco 3 ecm oll wea e210 aT O10 Ao1020 





+ vee + EL Coro 5020 5001) 
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where 


(= +1)! 9 
k, Qo1 G29 3 k even 
2 

Ee (0100015020) = 
+3) 

Dae (k-1)/2 

k Te Ok eae ! 9 k odd 

(F- 5)! 


Mewcan be shown Similarly. as stor the firstecol unne 


five t 
k/2] 
k-2m 
7 (k-mt1)! A106 

a » Con i a 
Bommei>s = )} (Column 3) 
AC1,4) = (541) aroaod 

: | + | + > : } 
moze) = gett yt2) THe) Bateag 

: * ap IL (Cg ieee i Gg fees 3 4 Cae Cane 2» 
mes, )) = yet ey rg Cge3) 19007 + a Samet @as 

mn V+] ae 7+) 4 (3 4+1)o 543) 2 : 
h(4,j) = Qe Di eh i gdod - aL 1 aie G15 Oe 


She 2k +2) 2 
uF : GoGo 


Zri 





[k/2] 


e Nee : 
mck, 3) = } (k-mt+j)! penal! O01 jon, 
Maco! Geeon nn 


3 Si (ORS) 


m=0 


The unit pulse response A is summarized for easy reference 


mielable 5.2. 


weaeeoteeG: All Absolute Entries of the 
Witter iwice Response 


The summation of all entries of <A> is performed 


following Cauchy's third method, which is outlined in 


Appendix K. 


The sum of all entries along the gee diagonal, where 


jemi, ...,l,... , using the notation of Appendix K and 


| 


D, = 1 (I 
Meee | 075! + |ao1|) Cine 
Do = Claro] + |aor])* + aao Cie 
Ds = Claro} + Joorl)? + 2a20laro| + Jaor|) lee 
ee laic| + laorld* + Sooo cvourc |e? ele ome) a 20° (I. 


T=2 


+ 


em= (!o10| lee we + (I-l)a2o€lario| + |aor]) 


seu) 


3b) 


3c) 


3d) 


3e) 


~, = = Z 
+ 5a PED? (lore + leans Tae + 
(ib S30) 
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It 1S proven in Appendix K that 


which can be written using eq (1.3) as 


+ 


aL 


+ 


Gra i, | + loo1|) 


Claro | a laoi|)? 


Claio| + laoi|)? 


© 
n=0 


n=0 
nev 
n=0 


Caco 


aGep 


n€n-1) ( n-2 


Q20) 


n(n-1)(n-2) n-3 


an Ca20) 


ele 





But this 1S equal to 


s = [1+ Clore] + fool) geo 


mrad g |. + area as 


O20 
3 q? 
+ Cloro] + Joaoil) ae cs | 
| 7 cae)" | Cia) 


te) 


where the | | designate a weighted infinite derivative 
SoemalOr. 
It 1s known from infinite series theory that the 


geometric series 


uf 


converges to =—_~—-- 
& Ll - Qo9 


af and Only 2f Vico gi <1) thus 
eq (1.4) can be rewritten, after applying the derivative 


SeeratoOr, as 


5 = 1 JY f doses teeth) © 
CL -a20) 7 


l - Qo 


om 
— 


= 
Ss 
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which converges to 


S$ = : ye ee 
~ (Cl = a20) 
122 ( Oo © ior ) 
1 - doo) 
merand Only if 
loio| + Jaorl < [1 - ao] Ci. 5) 
or equivalently, since laoq| <1, 
Q2 0 < 1l- Claro| + loo. |) 
Ge. Ge 


3. Improved Stability Conditions for Third-Order 
Uncoupled Filter in Two-Dimensions 


a 
Where Qe9 < 7 — 


It has been shown in eq (5.21) that the sum of all 


meme terms in the first column of Table 5.2 is: 


y In(n,o)| = » |(/=are )™ a a Ge) 
a 2¥-QA20 
n,=0 n,=0 


Si): 


ror 





2 
M10 
C20 mi 5) 
ipaere exists an angle 
O 
(watch oos as 


2¥Y-Q2 0 


for which 


Sg ( Q10 ) = ofl) (cos 6) 
] Nl BES SOL | 


But Since it can be shown [56]that: 


syalsel [Gel 22 seas!) 


Cr) ¥ 
One Weoc 1G) = =e 


it is possible to rewrite eq (1.6) in the following form 


epee 
q inGuoD | = » (=e) sintiny st) Oh Cr. 
ni=0 Tae) 


For this sum there exists no closed form solution. 
Following amethod by U. Apple [56], eq (1.7) can be viewed 
as the samples of a damped sine wave with continuous 


emepmment €etl)@ , starting at = -~l for which the 


E 
\ oe ES 


Bic 


ftest sample, mn, = € = -1 is, in this case, zero. This 
fact allows the determination of an upper bound by the 
equivalent integral, multiplied by a correction factor C 


as proven in [56]. 


& 
) Irmsol ses shy fc eam 
n,=0 =a) 
- |sin{(e +1)6]| de Cigse)) 


The argument of the integral has zero-crossings equally 


seagecd, e.£f£., 


i 5 -l<e< (itl): --1 were WSO, 25 oe 


ol 


The integral can thus be written as the integral over one 

segment (e.g., the first one) multiplied by the amplitude 
ws 

factor (vV-a29 + F and by summing over all i in the 


following manner: 
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oo 


€ 
| Caren.) jsanle +1)6lde 


alk 
TI 
go = . 
E i¢ 
Ss CY¥-ao9 ) sinfe +1)6 ae | > (¥-Qo20) 


al 1=0 


fe 0) 


a 
c (Gils V— O95 i) ) 


i; i 
eae. 
og (0 Fete Ce? 1 mbont ya) he)! 


tH 





1=0 


TT 
OC) + v~> O96 6 ) 


o}> 


Y-a29 (67 + [ln€¥-a20)]*) (1 - Yano 


moneda Only if Helm less ak : 


Thus an upper bound for 7 Cn OP es 
n,=0 


oo 


Tt 
6e(1 + Faro 6 ) 


q h(n, ,0)| < Ca ii eSeminar 
Far [67 + (In Fag )*) [1 - Fan FF) 
n ,=0 (1.9) 


where C is computed in Ref [56]to be C=1.1l1. 
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Theorem I: 


: q+1 
aes ) rca 5001 > » ain sab) 
n,; =0 n, =0 


mer each and every j. 
Proof: See Appendix L. 


The sum of all absolute entries in the unit sample 


response: 


— ‘a » [h(n ,n2) | 


n2=0 n,=0 


is bounded by 


10 yrs ae B 
1009) 
S » loo i | Dy rina 5001} 
n2=0 ni=0 
which by eq (I. 3) is bounded by 

ee T 

i +1 
N2 8 C1lt+v¥—ao9 : ) SE 
=an Y-O29 [07 + (In V-a20 )*] (1 - V-a20 FB) 


n2=0 
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which converges if 


aD] 3 


Sin 6° ¥-0120 (62 + In? ¥-Qoo0 J[1l- V-a20 ] 
Jaor] < 7 Cr 
es OGM Stara 6) 





0 
where Mec core cos: (==) 


2V¥-Q20 


Eeeeouh OF THE ABSOLUTE ENTRIES 

Peete UNTT SAMPLE RESPONSE 

Smee -OURTH ORDER (UNCOUPLED) FILTER 

The summation of the absolute entries of the unit sample 
response corresponding to the fourth order uncoupled transfer 
function follows Cauchy's third method, as outlined in 
Appendix K, and is Similar to the one performed in Section A 


of this appendix. 


Beewas shown in Section B.2.b of Chapter V, that 


©o rere) [n,/2]} ni-2i 7 [n2/2] oe 
6° a sit 20 a (nytn2-mei) !001 Cio 2 
q Ga aa Ven m! (n2-2m)! 
Mi-0 no-0 1=0 m=0 
Cina: 


1s equal to 
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[n;/2] Tee iio 2 


n2 my ~29 : as 
gC = la | lo | Q10 (n,tno-m1) mo 2 
ae a 1!{n,-21)! m! (no-2m)!} 2 
91 
1=0 m=0 


ni=0 n2=0 (Ei) 


mor Conditions stated in eqs (5.54) and (5.55). For con- 
venience, the unit sample response <A> is listed in 
Tables I.1 and 1.2, in which the absolute signs around each 
unit sample response entry, i.e., |h(n)| are replaced by 
absolute signs around 19 and around a9, in accordance 
ween eq (I[.11). 

The sum of all entries along the ie diagonal, l.e., 


after applying the binomial theorem becomes: 


i? 
we 1 
D, = Claro | uw |ao1|) 


Do = Cloro| + Jaor})* + azo + doe 
D3 = Claro| + Jaorl)? + 2€a20 + ao2)C]aio| + |aor]) 
Dy = Claro | i loor|)* + 3Cag9 t co2) Claro | |aoi|)? + (azo + do 2)* 
Ds = Claro | i |oor |)? + Ua, + do2) Claro | ‘ |ctor |)° a 
+:.3€a29 + O92)? Claro] + [aor] 
Dg = Claro tlaorj2® + Sazo + ao2)(Joro] + faor|?” 


t6Cao5 ao2)*Claro| + lena | + (29 + O92)? 


SO 


c 
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Fp ab Play eet ae : 
— 6 G\ on , fron 1p OI Cy ray 
plant eri PTE ea 
p z p Lk an e 
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piewsum Of all unit Sample response eMirie ce 


eS, 
a 


1=0 


is identified as 


Ss: = CQoo 2 re as 


oO 


Tee See 


+ Claro] + ioopa d 1h (eles) aE one a 


n=0 


ae Claro | 5 [ein |e y 


n=0 


But this is equal to 
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n(n=1) 


Zs 


(doo + Qo2) 


n-2 





But this is equal to 


: d 
S = ji + Claoi| + (Gl cles ar) 


a? 


dCaso 1 re 


=f 


2 
Stevo eo rl? 


‘ | 


! iy Caso x O55) 4 Ci le) 
n=0 
where | designates, Similarly as in the previous 


section, a weighted derivative operator. It is well known 


from geometric series theory that 


» C55 1 ae = eS 5 eee Cis) 
= Caso i Oo 2) 
n=0 
mmmemc Only if |ao5 * aoo| < 1. Considering this fact (1.10) 


and after applying the weighted derivative operator, eq (1.12) 


becomes 


BE i » } loro] * [aor] ? 
L = Caso + Ooo) d= Capa Qo 2) 


which converges to 


1 - Ca29 + doo) 1-{ SL ONGAN Pa 


1 - Cao6 u Oo 2) 


jmeeand only if 


loro] + |aorl <|2 - azo + aoa) | 
PMc thas can be rewritten since laanle: lage | <l as 
O20 + Go2 <1 - [aio] - [oor] 
Qe ud. 


See OUPLED FILTER: SECOND ORDER COUPLED (BILINEAR) FILTER 


Mnewcum Of absolute entries of A in Go lume Osis): 


Semverges to 


a : : 
Se Th eee ef andvonly 1f 7 ayo) < 1 


Bie 


Similarly , 


ee) 


= Ny n,-1l 
5} a | ICnitl ay Ooeen Marae o11| 


n,=0 


which can be rewritten, using derivative operator notation, 








as 
d n+l n 
= ) daze Caro “oor * oro a1) | 
n,=0 
cae n 
= ldare Coro%01 * Gia? » 10. 
n,=0 
| Se O10 ! Oe M10 aca i 
omer, = 2; 
wri +2 n 2 ni-l 
P i » (n din ) Barely + NiOj0° onl le) 
2 
oe nj(nj,-1) ~ n,-2 
7 g Gio O11 
a? 1. ni#2 ral en)...2 
7 os FZ O10 Qo1 t+ Ao Qo1011 + F O10 O11 
dour 9 
n,=0 
1 O10%01 + O11 2 | 
SS |) a EE li a © 
>|; - =o) IE Ci can | 


gil 


Her Ti> = 3: 


oO 
2 
31 0 SECS a Meera we O° Cloner 


2 yy | oat LL) nt 2) Ci to) ome, Se ay Fay tO eet 
= enn O 


n,=0 


4 (n,-1)n,(n,t1) ghtin4 4: (n1-2) (n}-1) (n)) ce 


DG 10 &91411 on 11 


3 


il 10a Bare? 
> 1001 et aa 
—_ 1 = Q19 1 - Qi0 


Following the same computational scheme, Ss becomes 


1 ioe oo J | 
Se > ——______ poset RL Dhol! ners Toe 
7 | 1 - Gio ! 1 - G10 =n | 
Since loro | < 1 fee. eo ) Secanebe 


written as follows 


i 
1 QioGo1i F 11 | 
S e > a ea en + a 
wo Lt = )6 } We Se a 
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APPENDIX K 


SUMMATION OF DOUBLE SERIES CAUCHY) 


Cauchy, in his note VIII of Analyse Algébrique and 
88 of the Résumés Analytiques [54], has defined for the 
two-variable case several methods of forming the sum of a 
double series, of which three are of importance for the 


derivation of stability conditions in two -dimensions. 
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First way: Sum to n terms each of the series formed by 
taking all terms in the first m rows of (K.1), where 


the sum of the first,m'" row is T sat 


Lyn PeSspecrlLvyeiy. 


ms 


Define 


Now, supposing each horizontal series to converge to 


ceca a and YT, to be a convergent series, where 


St os T = lim lim  §! 
m m,n 
moo N->0o 


m=0 
Second way: Sum to m terms each of the series in the 
j j U ee @ U e 
first n columns, and let the sums be iene ear 
Define 
ks = U ier aU 
m,n JIS n,m 


as a finite sum. When the vertical series converge to 
Ui,---5,U, respectively and )U, to be a convergent 


series, then 


- 
" 
= 
2 
it 


Lim Lim sit } 
Meso 
T->0o m->0o 
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Third way: Sum the terms Which lie in the successive 


diagonal lines of the array, namely AA',...,KK' and let 


the sums be See) respectively, i.e., 
Deo = paren 
= 28 
Ds Wyo * Wi 
= ay + 18 
Dat “ign Daal el 
Define a = Do et S.. + dD, as Velour tices um: Supposing 


Ps to be convergent, then 


Si = D 
Th 


n=2 


Cauchy's Theorem: Absolute Convergency of a Double Series 
————————— ee eeicy OF a Vouble Series 


Theorem: If all the horizontal series of 


{1 unm) |} 


be convergent, then 


Sst . : P 
: The horizontal series of 


= u(€n,m) 


ae 


al 


Sly 


are all absolutely convergent, and the sum of their sums 
momintinity Gonverges to a definite finite S. 
d, 


2 All vertical series converge absolutely and the 


sum of their sums to infinity converges to S. 
p40) The diagonal series is absolutely convergent 
and converges to 5S. 


Similar conclusions follow, if all the vertical or if 


the diagonal series of 


{| u(n,m) | 


be convergent. 


[Proof in Analyse Algébrique] 


ienean be concluded if one condition is found for which 
the entries of a double series converge absolutely, then 
the double series converges absolutely in every way the 


mathematical operation is performed. 
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AEE EN wie 1 


PeRO@E StHeOkREM APPENDIX I 


For 
ES Oyo 4 = Gg122 = O202 
let 
co co 
ez) = » » Gate An ee 
n,=0 n,=0 
Also 
il 
H(Z 152.) = MMM AG 
A 
giz.) (2 - Gey) 
; Notl 
2 iL n2 —-M2 
YY ieee 
ere 21 6 = C2041 
n,=0 


Bageting eq (L.2) and (L.2) results in 


yy q h(n, ,n.) 7 ae yy 3 hin, 0272] 


i, 0) gee eS agate 


Sik, 


Notl 


Glee!) 


(lie 2) 


Qo1 22 





Taking the absolute values of the corresponding coefficients 


and applying the triangle inequality leads for 21 = Zz =1 


iO 

notl ne 
y [h(ni yn2) | Jy y [h(n ,0) | | oo 1 | 
n,=0 n2=0 n,=0 no= 
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APPEND 


ROUGCHE op lUEOREMUAND RE ROGro Cr Char LE Rava 


The proofs in this appendix and the following conse- 
quence of Rouche's Theorem are made by applying methods of 


eerS)) 


Pee CONSEQUENCE OF ROUCHE'S THEOREM 


Mameorem Ml: Let f(z) =a + ayzt... F az fOr) i) oes 
be a polynomial with complex coefficients. Define 

)=bo + biz + + b Dat h ba = 8. ; 
mez) = Do ny? sans na , where , > a9 a. - an An. 3 
fmeeeOy...,n-l . Then f(z) is stable if and only if 


la < 1 and g(z) is stable. 


2 
Eee ROOF OF THEOREM 6.1 

MereeO(z) as defined in eq (6.15), then for eacen 
vector y = CY 19666 9 Vy? » where Y; is a complex number, 


qv (z) is defined as 
qv(z) = OC yz) 


Theorem M.2: The polynomial Q(z) is stable if and only if 


femal y , where max Ree pele ar(z) is stable. 


Bian 


3 





Proof, Theorem M2: (see also [73]) 

Q(z) is assumed to be stable. If a7 (2) is not stable 
for some Y , then q7( ez) =0 for the magnitude of 
miside the unit disc, i.e. |°z| soy dls Se delenpiete 
QC % 7) = Q OO ZY 1 ge009 °ZY QD) =n sae gs Uae | ozy. | Sole; “whiten sis ae 
@entradiction on the assumed stability of Q(z). 


Conversely, assume 


aoe) = QC¥125+++ YZ) 


is stable for all y , with maxly,| = 1. If 


QC °215++-5°2)) = O for Paani as 1h SUP tenet 


as-(r) = QCyi Ps... ¥yP) = 0 where pr = max| °z., | and 
aS 2. which is a contradiction. Gewds 
For z complex, and qji(z) = Q(z) , where Q is 


associated with CArs+++2A,) ; it follows that the unique 


leading coefficient of qi(z) is 


i.e., the sum of the leading coefficients of Q 


The polynomial q,(z) factors: 


Giez) = wCz=t1) *°s.(2 = ty) 


where Tioeseo ty are the complex roots. 
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leg (Ar y.++sA,) is stable, then Q(z) is, and hence it 


follows by the previous theorems that 








fees 1 LOT eelet= pile ae ie 
k 
poe | = hana) | = | w | | II t. 
peak 
I< 
where | II t. Pose |: 
ie 


Therefore |w| = | Vo A 
a 


See Or OF THEOREM 6.2 AND THEOREM 6.3 


This proof is very similar to the previous one. 


Let gz) = Q(P1 256-6 Py) [see eq 6.15] 
k 
55 
then q(1l) =1- AS (P1966 69 Py) 
Lad) 
k 
The ab 
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k 
jig » A;P. a 


then 


Ql ea 
mado) = O(0,...,0) = 1 


Hence by intermediate value theorem of one variable calculus 


there existS Zo 
ONES Ares oh 
suen that 
Gez = 0 


But IZ 0-7, | aecor + = Lowes were Pee cdo 
and this contradicts the hypothesis that (A1,..+5Ay,) and 


hence @ is stable. 
oracle be 


eee ROOF THEOREM 6.5 
The proof is similar to the previous ones. 


Let QC Zr 566+ 92) be the polynomial associated with 


Wyss es gA Define the polynomial Q'C 219-66 Zu) by 


ks 
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Q1CZ1 y+ 6+ Zu) = Q(061215-+-55y 2) 


teers NOt difficult to check that Q' is the polynomial 


associated with (Al,... Ay) Wa Ne AS are related to 


the A. as above. Of course, A.'S are real numbers if 
the Ps and A. are. Finally it is immediate that Q' 
[/mcrable if Q is, since |S, | =1. Therefore, 


(A1y-+-,A,) iometable of ) (Aja a 2S Stable, 


t 
k 


G2 end. 
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APPENDIX N 


PROObS “OF (Chaba Wal 


A. PROOF THEOREM 7 [See [73]). 

















Assume dg Zo) is unstable, 1.€., qxCZo) = 0 and 
fem) < 1. If Q(z) is stable, then 
6 
6 N 
DA a ees ao 
where > 1 for some i , then |zo| > = > 1 
i 


which is a contradiction. 
MmmOcrym lS Unstable 1t follows that there is a 


CO156- + 55y) where max |S. | Sl send “Cos. . -40,05— 0 


N 
IE ag(z) is stable, i.e. 





ag z) = Q & a 2 x) 


where r = max 6. | 
then agr) = Q 


Peer < t which is a contradiction. 


Pee ROOF OF COROLLARY 7 
Assume as in the previous proof that ds (Zo) 1s unstable, 


ae. dg Zo) =i. ston | Zo | <i 
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Meme rs stable, then 


Q(Z9515-+- 2054) = Q 
Where 
[Zo5.| > 1 for some i 


or 


ize = Sil 
| i 


Mien 2s a Contradiction. 


iio O(27) as unstable there exists a (6d ,...,6_) 


N 
such that max |6,| < 1 
i 


© 


et (515--. 563) = { 


we: ee a 
where r= max [6,| 
Let dg: (2) = QC2615--+525))) 


Then Agr (r) = QCS1,.-- 55) = 0 


ieee < 1 which is a COnERAGgT eC uo. 


Bye 
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